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ABSTRACT 

The  purpose  of  this  study  is  to  compare  two  noise  canceling  system 
algorithms  when  the  number  of  inputs  to  the  system  is  small.  The  expected 
value  algorithm  is  an  estimate  of  the  optimum  solution  for  minimum  mean  square 
error  analysis.  The  LMS  (least  mean  square)  algorithm  was  developed^  Widrow 
(Widrow  et  al_. , 1975)Ho  avoid  the  matrix  inversion  inherent  in  the  optimum 
solution  given  by  the  Wiener-Hopf  equation.  The  comparison  between  the  two 
algorithms  is  measured  by  the  rate  of  convergence  of  the  mean  square  error 
and  the  minimum  error  after  five  hundred  iterations.  ^ 

A two  input  noise  canceling  system  is  investigated  in  the  first  part 
of  the  study  and  a three  input  system  is  studied  in  the  second  part.  A two 
input  noise  canceling  system  has  one  input  that  contains  a mixture  of  the 
signal  and  noise  and  a second  input  which  is  noise  alone.  The  noise  only 

I 

input  is  multiplied  by  a weight,  W,  and  subtracted  from  the  first  input  to 

I 

produce  an  estimate  of  the  signal.  The  three  input  system  adds  another  noise 
only  input  and  a corresponding  weight. 

The  expected  value  and  LMS  algorithms  attempt  to  determine  the 
Wiener-Hopf  weight  vector  through  different  methods.  The  expected  value  al- 
gorithm estimates  the  necessary  expected  values,  substitutes  them  Into  the 
Wiener-Hopf  equation,  and  performs  the  necessary  matrix  Inversion.  The  LMS 
algorithm  uses  gradient  techniques  and  the  method  of  steepest  descent  to 
avoid  the  complexity  of  the  matrix  inversion. 
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A computer  simulation  of  a dc  signal  in  Gaussian  noise  provided  the 
comparison  between  the  two  algorithms.  The  tradeoff  In  terms  of  performance 
for  simplicity  in  the  LMS  algorithm  was  discovered  to  be  great.  The  ex- 
pected value  algorithm  consistently  produced  better  results  in  terms  of 
faster  rates  of  convergence  and  smaller  minimum  errors  after  five  hundred 
iterations.^ 

s \ 

The  LMS  algorithm  presented  problems  in  implementation  because  of 
two  parameter  values  that  had  to  be  chosen.  A value  for  the  Initial  weights 
had  to  be  chosen  since  the  weights  were  determined  using  only  previous  data. 
This  value  was  not  as  important  when  considering  such  a large  number  of  it- 
erations as  was  the  value  of  y,  the  term  that  controls  the  rate  of  converg- 
ence. This  term  must  be  within  a specific  range  in  order  to  guarantee  con- 
# 

vergence  of  the  algorithm.  This  range  is  small  when  only  considering  five 
hundred  Iterations.  Large  errors  produced  by  the  algorithm  were  also  a prob- 
lem. 

In  conclusion,  the  expected  value  algorithm  was  judged  to  be  the 
better  of  the  two.  The  differences  in  performance  were  great  especially 
when  the  slgnal-to-noise  ratio  or  the  correlation  between  the  noises  was 
small.  If  a "correct"  value  of  u is  chosen,  the  results  may  be  similar,  but 
if  not,  the  results  can  be  extremely  different  from  those  of  the  expected 
value  algorithm.  Thus  for  a small  number  of  inputs,  where  the  matrix  to  be 
Inverted  Is  of  small  dimension.  It  Is  better  in  terms  of  performance  to  use 
the  approximation  to  the  optimum  solution  given  by  the  expected  value  algor- 
ithm. 
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CHAPTER  I 
INTRODUCTION 

An  important  concept  in  signal  detection  theory  is  the  noise  cancel- 
ing system  which  is  used  to  eliminate  noise  from  a contaminated  signal.  The 
noise  canceling  system  consists  of  one  input  containing  the  contaminated 
signal  and  a varying  number  of  noise  alone  inputs.  Several  algorithms  have 
been  developed  in  order  to  decrease  the  number  of  computations  and  thus  in- 
crease the  speed  of  convergence  to  the  minimum  error. 

Early  work  in  least  squares  estimation  was  done  by  Wiener  in  1942 
when  he  recognized  that  the  solution  of  the  Wiener-Hopf  equation  determined 
the  optimum  estimate  for  stochastic  processes  (Kailath,  1974).  This  was  later 
extended  to  stationary  and  nonstationary  processes.  Kalman  also  contributed 
to  the  work  by  developing  a filter  which  modeled  the  covariance  of  the  sig- 
nal process  (Kailath,  1974).  The  LMS  (least  mean  square)  adaptive  algorithm 
was  developed  in  1959  by  Widrow  and  Hoff  at  Stanford  University  (Wldrow  et  al . , 
1975).  Work  on  adaptive  systems  continued  in  the  1960's  with  major  contribu- 
tions by  Lucky  at  Bell  Laboratories  (Wldrow,  et  al_. , 1975).  This  research 
compares  the  LMS  algorithm  with  an  estimate  to  the  optimum  solution  for  two 
and  three  Input  noise  canceling  systems. 

A general  n-input  noise  canceling  system  consists  of  one  input  with 
both  signal  and  noise  and  (n-1)  noise  only  Inputs.  Each  noise  only  input  is 


(2) 


multiplied  by  a weight  and  subtracted  from  the  first  input  in  order  to  ob- 
tain an  estimate  of  the  signal.  The  system  is  shown  in  Figure  1.1. 
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Figure  1.1.  An  n-input  Noise  Canceling  System. 


The  weights  are  chosen  so  that  the  signal  estimate,  s,  is  a linear 
least  mean  square  estimate.  That  is,  it  is  desirable  to  minimize 

EO2]  3 EC(s-s)2]  (1.1) 

through  a linear  combination  of  the  inputs.  The  optimum  solution  is  given 
by  the  Wiener-Hcpf  equation: 

W*  * R~\  (1.2) 

where  W*  is  the  optimum  weight  vector: 


(1.3) 


2 


X.  is  the  input  vector: 


4 


(1.4) 


R is  the  input  correlation  matrix: 


R = E[X  XT]  - 


EC^2x2^ 

E[x3x2] 


• £[x2xn] 


(1.5) 


[_E[xnx2]  ...  ECxnxn]J 

and  £ Is  the  cross  correlation  vector  between  the  desired  response,  d and  the 
input  vector  X_: 

(1.6) 

The  minimum  error  can  be  found  by  substituting  this  value  of  W_  into 
the  equation  for  the  error: 

*T 

s • Xi  - W X (1.7) 

and  thus 

min  E[(s-s)2]  * E[(x1  - W*  X-s)2]  (1.8) 


P - E[dX] 


E[dx2] 

E[dx3] 


E[dxnJ 


The  expected  value  and  IMS  (least  mean  square)  algorithms  attempt  to 
determine  the  Wiener-Hopf  weight  vector  through  different  methods.  The  ex- 
pected value  algorithm  estimates  the  expected  values  of  the  R and  £ matrices 
and  substitutes  these  into  the  Wiener-Hopf  equation.  Thus  the  weight  vector 
is  determined  by  actually  taking  the  matrix  Inverse.  The  LMS  algorithm 
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discussed  by  Widrow  et  al_. , in  "Adaptive  Noise  Cancelling:  Principles  and 
Applications"  (Proceedings  of  the  IEEE,  December,  1975)  uses  gradient  tech- 
niques and  the  method  of  steepest  descent  in  order  to  avoid  the  complexity 
of  the  matrix  Inverse  Inherent  In  the  optimum  weight  vector.  The  two  al- 
gorithms are  simulated  on  a POP-11/45  computer  at  Duke  University.  A plot 
of  the  mean  square  error  as  a function  of  the  number  of  iterations  provides 
the  main  basis  for  the  comparison. 

Chapter  II  derives  the  algorithms  for  the  two  input  noise  canceling 
system.  In  this  case,  the  R.  matrix  is  a 1x1  matrix.  The  optimum  weight 
vector  and  the  theoretical  minimum  mean  square  error  are  also  derived.  Chap- 
ter III  presents  the  computer  simulation  and  the  results. 

Chapter  IV  extends  both  algorithms  to  the  three  input  system.  The 
expected  value  algorithm  is  again  the  estimate  to  the  optimum  weight  vector 
formed  by  estimating  the  necessary  expected  values.  The  LMS  algorithm  is 
Widrow' s algorithm  where  the  one  element  weight  vector  of  the  two  input  sys- 
tem has  been  replaced  by  a two  element  weight  vector.  Chapter  V presents 
the  results  of  the  three  input  noise  canceling  system  computer  simulation. 

The  final  chapter,  Chapter  VI,  summarizes  the  results.  Directions 
for  possible  future  research  are  discussed. 
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CHAPTER  II 

DERIVATION  OF  THE  ALGORITHMS 


This  research  explores  the  problem  of  minimizing  the  mean  square  error 
of  a noise  canceling  system.  A two  input  noise  canceling  system  is  shown  in 
Figure  2.1. 


x1  a s + n1 


*2  3 n 2 


Figure  2.1.  Two  Input  Noise  Canceling  System. 


The  signal,  s,  is  uncorrelated  with  the  noises,  n^  and  n^,  which  are  correlat- 
ed with  each  other.  The  weight,  W,  attempts  to  transform  the  noise  n^  so  that 
It  can  be  subtracted  from  the  x^  input  to  give  an  estimate  of  the  signal. 

This  estimate  is  defined  as 


s ■ Xi  - WXg 


(2.1) 


Thus  the  error  is  given 
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The  mean  square  error  Is  defined  as  the  expected  value  of  the  squared  error: 

^ * E[e2]  - E[(s  - s)2]  (2.3) 

It  is  desired  to  derive  the  weight,  W,  In  such  a way  that  the  estimate 

of  the  signal,  s,  produced  is  a linear  least  mean  square  estimate.  Thus  the 

2 

estimate  attempts  to  minimize  the  mean  square  error  E[e  ] through  a linear 
combination  of  the  inputs.  The  optimum  solution  is  given  by  the  Wiener-Hopf 
equation.  This  problem  was  studied  by  Widrow  and  the  LMS  (least  mean  square) 
algorithm  was  developed.  This  algorithm  is  used  in  many  noise  canceling  sys- 
tems and  was  one  of  the  two  studied  in  this  research. 

1 . LMS  Algorithm 

The  LMS  algorithm  starts  with  the  basic  noise  canceling  structure 
shown  in  Figure  2.2. 


Figure  2.2.  Noise  Canceling  Structure  for  LMS  Algorithm 


The  LWS  algorithm  as  described  in  Widrow1 s "Adaptive  Noise  Cancelling: 

I 

Principles  and  Applications,"  attempts  to  adjust  the  weight  vector,  W,  to 
minimize  the  mean  square  error.  The  method  of  steepest  descent  is  used  as  a 
practical  way  to  find  close  approximate  solutions  to  the  Wiener-Hopf  equation: 


8 


I 


W*  * R_1P  (2.4) 

where  W*  is  the  optimum  weight  vector,  R is  the  input  correlation  matrix 

i-ECXjXj1]  (2.5) 

and  P_  Is  the  cross  correlation  vector  between  the  desired  response  and  the 
input,  or  X^,  vector 

P * E[d jXj]  (2.6) 

as  defined  in  Chapter  I for  the  jth  inputs. 

In  order  to  avoid  the  matrix  inversion  required  to  find  the  optimum 
weight  vector,  the  IMS  algorithm  sets  the  next  weight  vector,  , equal  to 
the  present  weight  vector  plus  another  term  proportional  to  the  negative 
gradient,  v.,  of  the  weight  vector.  Thus 

J 

■ Sl  - »2a  <2-7> 

where  y is  a parameter  that  controls  stability  and  the  rate  of  convergence. 
(For  a study  of  convergence  to  the  optimum  solution  see  Appendix  A).  The  in- 
stantaneous gradient  7j  is  also  estimated  by  assuming  that  the  squared  error 

2 2 
tj  Is  an  approximation  of  the  mean  square  error  E[ej  ].  Thus  the  LMS  algor- 
ithm uses  an  estimate  7.  of  the  true  gradient. 

J 

From  Figure  2.2,  the  estimated  signal,  s,  for  an  n input  system  is 

given  by 

s - x1  - WTX  (2.8) 

Substituting  (2.8)  into  the  value  of  the  error,  e,  given  by  (2.2)  yields 


« • - n 

Taking  the  partial  derivative  of  the  error  with  respect  to  W produces 

3C  _ y 

aff  * 


(2.9) 


(2.10) 
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(2.11) 


Ifi. 

aw 


The  partial  derivative  of  the  squared  error  with  respect  to  W is  given  by 

The  estimated  gradient  at  the  jth  Iteration,  v^. , Is  set  equal  to  the  partial 
derivative  of  the  squared  error  with  respect  to  W.  Thus 

2 

* -2^  (2.12) 

Finally  from  the  method  of  steepest  descent,  the  weight  vector  can  be  deter- 
mi ned : 

Jij+1  * 

-j+1  “ -j  + 2uej-j 
When  the  noise  canceling  system  has  only  two  inputs  as  shown  in  Fig- 
ure 2.3,  the  expression  for  the  weight  vector  reduces  to  a scalar  quantity: 

W,  • W-  +2ue,x,  (2.14) 

*j+l  J 2j 


(2.13) 


Figure  2.3.  Two  Input  Noise  Canceling  System. 

This  noise  canceling  system  can  also  be  written  in  terms  of  a feed- 
back system.  A general  two  Input  feedback  noise  canceling  system  Is  shown 
in  Figure  2.4.  The  input  d^  is  the  desired  response  and  cj  is  the  error. 


Figure  2.4.  General  Two  Input  Feedback  Noise  Canceling  System. 

Since  the  desired  input  is  actually  the  unknown  signal,  dj  is  often 
set  equal  to  x-j  (Widrow  et  al_. , 1975).  When  this  is  done,  the  system  is  not 
a feedback  system  in  the  true  sense  because  the  error  term  is  now  given  by 

ej  * dj  " slMSj 

ej  * xij  - s\ms. 

ej  * WLMSjX2j  ^2*15: 

Thus  the  weight,  Wj>+1 , is  given  in  terms  of  inputs  only  when  (2.14) 
is  used  to  define  the  weight. 

Since  the  noise  canceling  system  employing  the  LMS  algorithm  as  de- 
fined by  (2.14),  is  not  an  actual  feedback  system,  another  algorithm  was 
developed.  This  algorithm,  the  expected  value  algorithm,  also  attempts  to 
minimize  the  mean  square  error.  However,  a feedback  approach  was  not  taken. 
Instead  two  weights  and  a brute  force  method  of  solving  the  matrix  Inverse 
avoided  by  the  LMS  algorithm  were  used.  For  a two  input  system  this  matrix 
is  one  dimensional  so  the  Inverse  is  quite  simple. 


2.  Expected  Value  Algorithm 

The  expected  value  algorithm  uses  a ratio  of  estimates  to  the  ex- 
pected values  of  the  inputs  in  order  to  minimize  the  mean  square  error.  The 
weight,  W,  is  derived  for  the  expected  value  algorithm  using  the  two  Input 
noise  canceling  structure  shown  in  Figure  2.5. 

x,  * s + n,  ^ + x- — v 


Figure  2.5.  Two  Input  Noise  Canceling  System. 


For  this  system,  the  estimated  signal,  s,  is  given  by 
s * ax.j  + bx2 

The  orthogonality  principle  gives 


E[(s  - sjx^  » 0 
E[(s  - s)x2]  » 0 


(2.18 


Substituting  the  value  of  s given  by  (2.16)  into  (2.17)  and  (2.18)  produces 


a£[x-|X^]  + bECx^]  - E[sx]  ] » 0 

and  a£Cx-|X2]  + bECx2x2]  - E[sx2J  ■ 0 

Writing  (2.19)  and  (2.20)  In  matrix  form  gives 

fECxlXl]  ECxgX^I  fal  [ECsx^ 


E[xnx,]  E[x,x,]  b E[sx,] 


(2.19 

(2.20 


(2.21 
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Taking  the  matrix  inverse,  equation  (2.21)  can  be  solved  for  a and  b.  Thus 


and 


ECsx^Efxg  ] - ECx^jECsXg] 
1 E[x,Z]E[x/]  - [ECx^j]]2 


ECx1 2]E[sx2D  - E[x1x2]E[sx1] 
E[x12]ECx22]  -~[E[x1x233k 


(2.22) 


(2.23) 


By  hypothesis,  E[sn^]  = 0 and  E[sn2]  3 0.  Thus  E[sx2]  =»  0 and 
E[sx-|]  * E[s(s  + n1 )]  = E[s2]  + E[sn.,]  * E[s2] 


'1 


1- 


(2.24) 


Substituting  these  values  into  (2.22)  and  (2.23)  gives 

I? 


a * 


7 


E[s2]E[x22] 
T- 


e[x/]e[x24]  - CeCx1x233" 


(2.25) 


and 


E[s^]E[x1x2] 

ECx12]E[x22]  - [E[Xix2]]2 


(2.26) 


Let 


„ E£sfj 

E[x,2]E[x22]  - [E[x,x2]]2 

Substituting  the  value  of  a into  (2.25)  and  (2.26)  gives 


(2.27) 


a » aE[x/>  ] 


and 


(2.28) 

(2.29) 


b * -aECx-jXg] 

Looking  again  at  the  noise  canceling  system  diagram  and  replacing  the 
values  of  a and  b as  determined  in  (2.28)  and  (2.29)  yields 


x,  * s + n 


Figure  2.6.  Two  Input  Noise  Canceling  System  for  the  Expected 
Value  Algorithm. 


s + n 


Figure  2.7.  Two  Input  Noise  Canceling  System  for  the  Expected 
Value  Algorithm. 


The  term  a can  be  simplified  by  recognizing  several  relationships 
One  input  to  the  noise  canceling  system  is  x-j  where 

X1  * s + n,  ( 

Squaring  both  sides  of  (2.30)  gives 

x-j2  * s2  + n^2  + 2sn^  ( 


Since  E[sn.j]  • 0,  taking  the  expected  value  of  both  sides  of  (2.31)  gives 
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■ 


ECx,2]  * E[s2]  + ECn-,2]  (2.32) 

The  second  input  to  the  noise  canceling  system  is  x2  where 

x2  - n2  (2.33) 

2 2 2 7 

Since  x2  = n2  , taking  the  expected  value  of  both  sides  gives  E[x2  ]»E[n2  ]. 

Also  since  x^x2  3 sn2  + n^n2>  taking  the  expected  value  yields  Efx^]  * 

ECn^]. 


The  correlation  between  n-j  and  n2  is  defined  as  where 

E[n-jh2] 


012 


VeCn^jECn.,2] 


(2.34) 


Squaring  and  rearranging  (2.34)  produces 

[E[nin2]]2  = P122E[n12]E(n22] 


(2.35) 


Substituting  (2.35)  into  (2.27)  gives 


a 

Since  E[x22] 


E[s21  ] 3— 

ECn22 ] [E[s2]  + ECn/]]  - p^E^2] 

ECn22],  let 

0 3 aE[x22]  3 ot£[n22] 


(2.36) 


Thus 


its!l 


E[s2]  + E[n.2](l  - o 


12 


(2.37) 


The  system  of  Figure  2.7  can  be  reduced  again  to  produce  a multiplier 
of  one  in  the  top  branch  as  shown  in  Figure  2.8. 

The  term  s is  found  to  depend  on  the  apriorl  knowledge  in  terms  of 
the  slgnal-to-nolse  ratio  and  the  correlation  between  the  two  noise  terms. 
Oividing  the  numerator  and  denominator  of  (2.37)  by  E[n^  ] and  recognizing 
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Figure  2.8.  Two  Input  Noise  Canceling  System  for  the  Expected 
Value  Algorithm. 


that  the  signal-to-noise  ratio  (S/N)  is  the  term  E[s2]/E[n^]  yields 

$ 

N 


W+  1 - p 


(2.38) 


12 


Thus  the  weight,  W,  for  the  expected  value  algorithm  is  given  by 

E[x,xj] 


ECxJ] 


(2.39) 


where  the  term  s must  also  be  a multiplier  of  the  output  of  the  system. 

In  sumnary,  the  expected  value  algorithm  uses  a ratio  of  the  expected 
values  of  the  inputs  in  order  to  determine  the  weight,  W.  However,  since 
expected  values  are  mathematical  concepts,  they  are  estimated  for  the  com- 
puter simulation  In  terms  of  the  jth  sequence  of  inputs.  The  two  input 
noise  canceling  system  with  the  expected  value  algorithm  is  shown  in  Figure 


2.9. 

Both  the  LMS  and  expected  value  algorithms  have  only  one  weight 


which  is  multiplied  by  the  noise  only  input.  Therefore  an  Interesting  prob- 
lem is  to  find  the  optimum  weight  for  a linear  least  mean  square  estimate. 


Figure  2.9.  Expected  Value  Algorithm. 

Optimum  Weight  Vector 

The  noise  canceling  structure  is  shown  in  Figur?  2.10. 


"I 

'2  * n2 



Figure  2.10.  Two  Input  Noise  Canceling  System. 

The  error,  e,  for  a linear  least  mean  square  estimate  is  given  by 

e - E[(s  - s)2]  (2.40) 

where  s Is  defined  in  (2.1).  Substituting  the  value  of  s given  by  (2.1) 

Into  (2.40)  yields 

c - E[(Xl  - Wx2  - s)2]  (2.41) 

Taking  the  derivative  of  the  error  with  respect  to  the  weight  and  setting  the 


result  equal  to  zero  produces  the  optimum  weight.  Thus  taking  the  partial 
derivative  of  (2.41)  with  respect  to  W produces 

» E[(x1  - W x2  - s)2]  • 0 


and  thus 


8y  hypothesis 


ECx1x2]  - WE[x2Z]  - E[sx2]  - 0 


ECs^]  » E[sn2]  • E[sx2]  ■ 0 
Solving  for  W in  (2.42)  yields 

E[x,x?] 

Wnnt  - '-f-  (2.4 

opt  ECx22] 

Thus  the  optimum  weight  for;  a linear  least  mean  square  estimate  is 
given  by  (2.43)  and  the  noise  canceling  system  is  shown  in  Figure  2.11.  Thl 
weight  is  the  same  as  that  given  by  the  Wiener-Hopf  equation  for  a two  input 
noise  canceling  system. 


x,  * s + n. 


Figure  2.11.  Optimum  Weight  Vector  for  the  Two  Input 
Noise  Canceling  System. 
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Notice  that  the  optimum  weight  and  the  theoretical  weight  for  the  ex- 
pected value  algorithm  are  the  same.  The  expected  value  algorithm,  however, 
has  a multiplier  of  8 which  is  not  present  in  the  noise  canceling  structure 
of  Figure  2.11.  When  8*1,  the  expected  value  algorithm  would  be  the  same 
as  the  optimum  if  the  expected  values  could  be  computed  directly.  The  weight 
for  the  IMS  algorithm  will  converge  to  the  optimum  weight  as  the  number  of 
iterations  increase  (see  Appendix  A). 

Since  both  the  IMS  and  expected  value  algorithms  have  weights  which 
converge  to  the  same  value,  the  theoretical  minimum  mean  square  error  can  be 
found. 

4.  Theoretical  Minimum  Mean  Square  Error 

The  theoretical  minimum  mean  square  error  for  both  algorithms  is  the 
same  (since  the  weights  converge  to  the  value  given  by  the  optimum  weight) 
and  can  be  determined  from  the  optimum  weight  vector.  Using  the  noise  can- 

( 

celing  structure  of  Figure  2.10  and  the  optimum  weight  vector  given  by  (2.43), 
the  minimum  mean  square  error,  e^,  can  be  computed: 

i ‘ms  ' £C(s  - s)2] 


and 


If 


t*B’£C(xi ' luTT*2  ■ s>  ] 

ECn,n_3  ? 

eMMS  * E^nl  " -,“2,  n2^  ^ 

ELn2  J 


Y 


E[n1n23 

" ECn22] 


(2.44) 


(2.45) 


then  substituting  (2.45)  into  (2.44)  produces 


Expanding  (2.^)  yields 


*MMS  * 3 * ^yECn-jng]  + y EC^  ^ (2.4 

The  correlation  coefficient,  p12»  is  defined  in  (2.34)  and  can  be  re 
written  in  terms  of  ECn^].  That  is: 

ECn1n2^  * P1 2 J E^nl  2^E^n2  ^ 


Substituting  (2.48)  into  (2.45)  yields 


T, 


E[n<i  ] 


Y * P 


12, 


EC",2] 


Using  (2.48)  and  (2.49)  the  minimum  error  can  be  simplified: 


(2.4 


(2.4 


‘m  * E[ni  ] 
and  thus 


E[n/] 


*121  ECng7]  j 


°uJE[nlZlE[n22]  + 


Zil 


EC"!  ] 


’12  ECn^ 


ECn/] 


eMMS  * ^"l^1  " °12  ^ 


(2.51 


The  theoretical  minimum  mean  square  error  for  the  noise  canceling 
system  of  Figure  2.10  is  given  by  (2.50).  The  minimum  error  depends  only  on 
the  power  of  the  noise  mixed  with  the  signal  and  the  correlation  between  the 
two  noises.  As  the  correlation  decreases  to  zero,  the  error  approaches  its 
largest  value.  When  the  correlation  actually  equals  zero,  this  is  the  worst 
case  in  terms  of  error.  Thus  the  error  would  be  the  same  If  no  noise  cancel 
Ing  system  had  been  used.  When  the  correlation  equals  one,  the  error  equals 
zero  and  the  best  case  is  realized.  Thus  as  the  correlation  Increases  and 
the  noise  power  decreases  (implying  Increasing  signal-to-noise  ratio),  the 
theoretical  minimum  mean  square  error  also  decreases. 
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For  the  expected  value  algorithm,  Figure  2.10  assumes  that  s * 1. 

This  is  a reasonable  assumption  because  8 depends  on  apriori  knowledge  of  the 
signal -to-noise  ratio  and  the  correlation  coefficient  p12*  IT  this  know- 
ledge is  known,  a better  estimate  can  be  obtained. 

Using  the  noise  canceling  structure  of  Figure  2.9  and  the  value  of 
8 from  (2.38),  the  theoretical  minimum  mean  square  error  for  e f 1 can  be 
computed.  From  Figure  2.9,  the  estimated  signal,  s,  is  given  by: 


E[x,x2] 

s * 8[X1  • Trrrr  x2] 


E[x  ‘] 


(2.51) 


Substituting  (2.51)  into  (2.40)  gives 


"MMS 


E[x,x?]  2 

E[(8x,  - 2 — 0X2  ' 

1 e[x22]  2 


(2.52) 


Using  (2.43)  and  (2.48),  equation  (2.52)  can  be  simplified: 

eMMS  * ECni2](1  - °i22>a2  + ECs2](8  - I)2 
If  8 3 1,  then  (2.53)  reduces  to  (2.50). 


(2.53) 


5.  Summary  of  the  Algorithms 

The  algorithms  were  derived  in  a theoretical  manner,  but  are  summarized 
here  in  the  forms  which  were  implemented  for  the  computer  simulation.  The  in- 
puts, weights,  and  estimates  are  given  In  terms  of  the  jth  iteration. 

A.  IMS  Algorithm 

The  weight  vector  for  the  IMS  algorithm  was  defined  in  (2.13)  In  terms 
of  the  error: 

(2.13) 

However,  since  the  signal,  s,  is  not  generally  known.  It  is  desirable  to  re- 
place the  error  cj  in  (2.13)  by  §j.  Then  the  weight  vectors  can  be  shown  to 
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converge  to  the  same  value.  Thus  for  the  jth  Iteration 

CJ  * ' SJ 

Replacing  by  i.  in  (2.13)  gives 


*i,+,  ‘ * 2“‘A  * 2“SA 


and  in  general 


(2.54) 


(2.55) 


j-i  j-i 

s*  ■ So  * 2» E *A  *2»I> 

J i*0  i*0 


Taking  the  expected  value  of  both  sides  of  (2.56)  produces 

.1-1  M 

E[W*  3 - Wq  + 2u  Y]  E[eiX.]  + 2u  V E[X.Si] 
j 1-0  1-0 

By  hypothesis  E[X_  s]  * 0.  Therefore 


(2.56) 


(2.57) 


J-i 

k ] ■ So + 2“  Z 

J 1-0 


(2.58) 


Therefore  replacing  ej  by  either  §j  or  3 s^.  - s^.  produces  the  same  weight 


solution. 


The  weight  vector  used  throughout  the  remainder  of  this  paper  re- 


places cj  by  Sj  and  for  the  two  input  system  is  given  by 


Wj+1  " MJ  + 2wV2, 


(2.59) 


The  weight  vector  can  also  be  written  in  terms  of  the  data  inputs  alone  by 


substituting  for  the  value  of  Sj  from  (2.1).  Thus 


Wj+1  " V1  “ ZuX2j  } + 2wXljX2j 


(2.60) 


The  weight  vector,  Wj+^ , is  computed  during  the  computer  simulation 
after  determining  Sy^.  This  value  is  stored  in  a vector  so  that  it  can  be 


used  to  compute  sLMS  for  the  next  iteration.  Figure  2.12  shows  the  IMS  algor- 
ithm for  the  jth  iteration  of  the  computer  simulation. 


B.  Expected  Value  Algorithm 

When  the  expected  value  algorithm  was  implemented  on  the  PDP-11/45 
computer,  the  expected  values  were  estimated  using  an  iterative  approxima- 
tion. For  example: 

j 

Zv2i 

ECxl  .x2  ] " (2.61) 

J J 


The  term  a was  set  equal  to  one  for  reasons  explained  later.  Thus  for  the  jth 
iteration,  the  expected  value  algorithm  is  implemented  as  shown  in  Figure  2.13. 

The  expected  value  algorithm  can  be  written  in  a form  which  reveals 
the  way  in  which  the  weight  is  updated: 


j 


I>i,V 
1-1  1 1 


WM+V*i 

— — -v- 


Wj-1  + x2 


j 


(2.62) 
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This  form  also  indicates  the  estimation  of  the  expected  values  needed  for 
the  optimum  weight  vector  defined  in  (2.43)'. 


C.  Compari son  of  the  A1 qori thms 

The  optimum  weight  vector  as  given  in  Figure  2.11  will  produce  the 
best  results  In  terms  of  mean  square  error  estimation  for  the  two  input  sys- 
tem. However,  since  the  expected  values  are  not  known  quantities,  they  must 
be  estimated.  The  difference  in  performance  between  the  expected  value  and 
IMS  algorithms  can  be  traced  to  the  weights.  The  expected  value  algorithm 
makes  two  estimates  (one  for  the  numerator  and  another  for  the  denominator) 
before  dividing  the  two  to  obtain  an  estimate  of  the  weight: 


Wj-1  + xl/2, 

— — v- 


Wj-1  + x2 


j 


(2.63) 


The  jth  Inputs  are  used  to  update  the  (j-l)st  estimates  before  the  division 


is  made. 

The  IMS  algorithm  however,  makes  one  estimate  for  the  weight: 
“LMSj  ■ “LHSj.,  * 
or  in  terms  of  the  inputs: 
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(2.64) 


W,  Mc  3 Wi  Mc  (1  - 2ux-  + 2ux,  x-  (2.65) 

LMSj  LMSj-l  2j-l  'j-l  2j-l 

2 

Although  the  same  terms,  x^x^  and  x^  , are  being  calculated,  the  two 
estimate  approach  of  the  expected  value  algorithm  generally  will  produce  a 
better  estimate  in  fewer  iterations.  The  expected  value  algorithm  also  uses 
the  jth  Inputs  to  compute  the  jth  weight,  but  the  LMS  algorithm  only  uses 
the  data  through  the  (j-l)st  inputs.  This  is  the  reason  the  IMS  algorithm 
can  be  written  as  a feedback  system  and  the  expected  value  algorithm  can  not. 


CHAPTER  III 

COMPUTER  SIMULATION  OF  TWO  INPUT  NOISE  CANCELING  SYSTEM 

The  two  algorithms  described  and  developed  in  Chapter  II  were  simu- 
lated on  a PDP-11/45  computer  at  Duke  University.  In  order  to  make  the  al- 
gorithms more  comparable,  the  noise  canceling  structure  of  Figure  2.10  was 
used  for  the  simulation.  The  multiplier  of  S for  the  expected  value  algor- 
ithm, as  shown  in  Figure  2.8,  was  set  equal  to  one.  If  a is  actually 
calculated,  a better  estimate  could  be  obtained.  However,  for  purposes  of 
comparison  it  was  desirable  to  use  the  same  amount  of  apriori  knowledge  for 
both  algorithms.  The  comparison  was  based  on  a graph  with  the  mean  square 
error  plotted  as  a function  of  the  number  of  iterations. 

1 . Problem  Description 

Values  for  the  inputs  to  the  noise  canceling  system  of  Figure  2.10 
had  to  be  chosen  before  beginning  the  computer  simulation.  A dc  signal  was 
used  as  the  signal  component  of  the  system.  The  value  of  the  signal  was  con- 
stant and  arbitrarily  chosen  to  be  10.0.  The  choice  of  a dc  signal  enabled 
a more  comprehensive  study  since  a frequency  dependent  signal,  such  as  a 
sine  wave,  would  have  had  the  effects  of  changes  In  phase  or  frequency  em- 
bedded In  the  results. 

The  noise  was  generated  using  the  computer's  random  number  generator. 
Uniformly  distributed  random  numbers  were  obtained  and  then  transformed  to 

(25) 
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zero  mean  Gaussian  random  numbers  by  employing  the  Box-Muller  equations. 

Each  pair  of  Gaussian  random  numbers  had  a given  variance  and  correlation  and 
served  as  the  noise  components  needed  for  the  inputs  to  the  noise  canceling 
system. 

The  algorithms  were  implemented  on  the  computer  using  several  iter- 
ative estimates.  The  expected  values  needed  for  the  expected  value  algorithm 
were  computed  using  the  approximation 

j 

SXliX2i 

EC*-,  x-  ] 3 : (3.1) 

J j J 

The  mean  square  error  was  also  approximated: 

j 


««Ji  - ■ “-j <3-2> 

The  initial  weight  vector,  W(l),  for  the  IMS  algorithm  was  set  to  a 
constant  value  of  5.0.  This  value  was  an  arbitrary  choice.  Since  the  mean 
square  error  was  computed  for  five  hundred  iterations,  the  effect  of  the 
initial  weight  vector  washed  out.  If  it  was  necessary  to  use  a smaller  num- 
ber of  iterations,  the  effect  of  the  initial  weight  vector  would  be  more 
important  and  thus  it  should  be  carefully  chosen. 

The  LMS  algorithm  also  required  values  for  the  parameter  u which  con- 
trols the  rate  of  convergence  of  the  algorithm.  In  order  for  the  two  input 
system  LMS  algorithm  to  converge,  the  value  of  u must  be  within  a given 
range  (see  Appendix  A)  defined  as 

|1  - 2uE[x22]|  < 1 (3.3) 

As  the  value  of  u decreases,  the  rate  of  convergence  also  decreases  since 
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|1  - 2u£[x2  ]|  approaches  1.  Three  values  of  y were  chosen  for  the  computer 
simulation.  These  were  values  of  .01,  .001,  and  .0001  which  correspond  to 
graphs  b,  c,  and  d,  respectively  on  each  page  of  figures.  For  convenience 
in  the  computer  simulation 

EC^2]  3 E[n22]  (3.4) 

Therefore  when  the  signal-to-noise  ratio  is  10.0,  all  three  y values  are 
within  the  range  for  convergence  of  the  algorithm.  When  the  signal -to-noise 
ratio  is  1.0,  y values  of  .001  and  .0001  are  within  the  range,  but  a y value 
of  .01  results  in  equality  in  equation  (3.3).  When  the  signal-to-noise 
ratio  is  .1,  only  the  y value  of  .0001  is  within  the  range  of  convergence. 
Equality  in  equation  (3.3)  results  for  a y value  of  .001  and  a y value  of 
.01  is  outside  the  range  for  convergence  of  the  algorithm. 

The  theoretical  minimum  mean  square  error  is  defined  in  (2.50)  and 

2 

depends  only  on  the  values  of  E[n^  ] and  p^2.  Thus  when  p^2  = 0,  the  worst 
error  is  produced.  Figure  3.1  shows  the  worst  case  for  the  signal-to-noise 
ratio  of  10.0.  When  p^2  3 1,  the  best  case  is  realized  and  the  noise  is 
completely  canceled.  Therefore  when  p^2  is  between  0 and  1 in  value,  some 
noise  canceling  is  achieved,  but  some  error  also  exists. 

With  the  noise  canceling  system  described,  the  method  of  comparison 
was  a graph  of  the  mean  square  error  as  a function  of  the  number  of  itera- 
tions. Each  curve  consisted  of  five  hundred  iterations.  Twenty-five  curves 
were  averaged  point  by  point  before  this  ensemble  average  was  plotted.  In 
this  way,  the  effect  of  a single  ensemble  (or  particular  sequence  of  random 
numbers)  was  less  and  the  averaged  curve  was  a smoother  curve.  The  curves 
can  be  compared  directly  or  the  difference  from  the  theoretical  minimum  mean 
square  error  after  five  hundred  iterations  can  be  used  as  a basis  of  compar- 
ison. 
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2.  Performance  Compari son 

The  graphs  of  the  mean  square  error  as  a function  of  the  number  of 
Iterations  follow.  They  are  arranged  by  pages  with  the  same  signal -to- 
noise  ratio  (S/N)  and  noise  correlation  (p.^)  on  each  page.  Each  page  con- 
tains a graph  using  the  expected  value  algorithm  and  three  graphs  of  the 
LMS  algorithm,  each  of  which  employs  a different  y value.  The  pages  are 
arranged  in  order  of  decreasing  signal-to-noise  ratio  and  decreasing  corre- 
lation coefficients.  All  of  the  graphs  were  plotted  using  the  same  scale  in 
order  to  compare  them  directly.  When  the  mean  square  error  was  too  large  to 
be  within  the  bounds  set  by  the  scale,  the  curve  was  plotted  as  a straight 
line  until  the  iteration  when  the  value  came  within  bounds.  This  does  not 
mean  that  the  mean  square  error  was  a constant  value  through  this  Iteration, 
only  that  it  was  too  large  for  the  desired  scale  (see  Appendix  C).  Table 
3.1  gives  representative  values  of  the  curves. 

The  difference  between  the  minimum  mean  square  errors  of  the  curves 


Is  explained  by  the  different  rates  of  convergence.  The  theoretical  minimum 
error  is  the  same  for  each  set  of  values  of  signal-to-noise  ratio  and  noise 
correlation.  Changes  in  the  amount  of  noise  correlation  have  little  effect 
on  the  rate  of  convergence  itself,  but  do  change  the  value  of  the  theoretical 
minimum  error. 

Figure  3.1  demonstrates  the  worst  case  (p^  * 0*0)  ^or  a signal-to- 
noise  ratio  of  10.0.  The  expected  value  algorithm  shows  rapid  convergence 
to  the  theoretical  minimum  error  of  10.0.  After  five  hundred  iterations, 
the  actual  error  (with  a value  of  11.22)  is  very  close  to  the  theoretical. 
This  Is  the  worst  case  for  this  signal-to-noise  ratio.  Any  larger  value  of 
noise  correlation  will  Improve  the  ability  of  the  system  to  eliminate  the 
noise.  The  LMS  algorithm  does  not  perform  as  well  as  the  expected  value 


29 


algorithm  because  the  rate  of  convergence  is  slower  and  the  minimum  error 
after  five  hundred  iterations  is  larger.  As  the  value  of  u is  decreased,  the 
rate  of  convergence  becomes  less.  However,  the  u value  of  .001  produces  a 
smaller  error  than  does  the  u value  of  .01.  Although  the  rate  of  convergence 
is  less,  the  minimum  error  is  smaller  due  to  fewer  fluctuations  in  the  weight 
vector. 

Figures  3.2  through  3.4  all  have  the  signal-to-noise  ratio  of  10.0. 
Each  figure  has  a different  value  of  noise  correlation.  For  each  figure,  the 
expected  value  algorithm  exhibits  the  best  performance  in  terms  of  a lower 
minimum  value  and  a faster  rate  of  convergence.  In  each  case  the  LMS  algor- 
ithm with  u equal  to  .01  compares  more  favorably  with  the  expected  value  al- 
gorithm in  terms  of  rate  of  convergence.  However,  with  y equal  to  .001,  the 
LMS  algorithm  compares  better  in  terms  of  minimum  error.  The  u equal  to 
.0001  case  of  the  LMS  algorithm  has  a slow  rate  of  convergence  and  large  ini- 
tial and  final  errors.  Thus  this  case  is  not  as  important  in  the  comparison 
between  the  two  algorithms. 

When  the  signal-to-noise  ratio  is  10.0,  small  changes  in  correlation 
produce  only  small  changes  in  the  curves.  As  the  correlation  is  reduced, 
the  rate  of  convergence  of  all  of  the  curves  decreases.  The  minimum  error 
after  five  hundred  iterations  becomes  larger.  The  maximum  error  for  the  ex- 
pected value  algorithm  never  exceeded  100,  but  the  LMS  algorithm  often  pro- 
duced large  errors.  As  the  correlation  decreased,  these  maximum  errors  in- 
creased. 

When  the  signal-to-noise  ratio  was  lowered  to  1.0,  the  comparison  is 
more  dramatic.  Figures  3.5  through  3.7  display  these  results.  Again  the  ex- 
pected value  algorithm  produces  the  best  results.  The  change  in  signal-to- 
noise  ratio  has  scarcely  affected  the  performance  of  this  algorithm.  However, 
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the  IMS  algorithm  has  been  affected  by  the  change  to  a great  degree.  The 
u equal  to  .01  case  of  the  IMS  algorithm  produces  a curve  that  is  entirely  off 
of  the  scale.  This  is  due  to  the  fact  that  the  value  of  u is  outside  of  the 
range  of  convergence.  The  u equal  to  .0001  curve  has  a very  slow  rate  of 
convergence  and  the  Initial  values  are  so  large  that  the  straight  line  effect 
is  observed.'  Thus  the  u equal  to  .001  curve  is  the  only  curve  for  the  LMS 
algorithm  that  can  be  effectively  compared  with  expected  value  algorithm.  The 
rate  of  convergence  is  much  slower  than  for  the  expected  value  algorithm  and 
the  minimum  error  is  much  larger.  The  LMS  algorithm  again  produces  large 
initial  errors. 

With  a smaller  signal-to-noise  ratio,  changes  in  noise  correlations 
are  more  Important.  Small  differences  in  correlation  have  a greater  effect 
on  the  rate  of  convergence  and  the  maximum  and  minimum  errors  than  when  the 
signal-to-noise  ratio  Is  larger.  If  the  correlation  is  constant  and  the 
signal-to-noise  ratio  is  lowered  or  if  the  signal-to-noise  ratio  is  constant 
and  the  correlation  decreases,  the  rate  of  convergence  is  slower  and  the 
maximum  and  minimum  errors  for  five  hundred  iterations  Increase.  This  was 
true  for  all  of  the  cases  tested  except  that  the  maximum  error  for  the  ex- 
pected value  algorithm  was  never  greater  than  100. 

When  the  signal-to-noise  ratio  was  lowered  to  .10,  the  expected  value 
algorithm  alone  produced  curves  with  values  small  enough  to  be  within  the  de- 
sired scale.  The  u values  of  .01  and  .001  are  too  large  for  the  algorithm  to 
converge.  However,  as  the  u valje  gets  smaller,  the  rate  of  convergence  be- 
comes too  slow  for  convergence  within  five  hundred  iterations.  Thus  the  u 
value  of  .0001  also  produces  a curve  with  values  too  large  for  the  scale. 

Small  changes  in  correlation  produce  noticeable  changes  in  the  expected  value 
algorithm  curves.  When  the  signal-to-noise  ratio  was  10.0,  a change  in 
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correlation  from  .99  to  .98  produced  barely  distinguishable  changes  in  the 
curves.  Nevertheless,  the  expected  value  algorithm  produces  good  results  in 
terms  of  a fast  rate  of  convergence  and  a minimum  error  that  approaches  the 
theoretical  minimum. 

3.  Summary  of  Two  Input  Results 

In  all  of  the  cases  tested,  the  expected  value  algorithm  produced  the 
superior  results.  The  expected  value  algorithm  consistently  had  a faster 
rate  of  convergence  and  a smaller  minimum  value.  It  also  had  two  fewer  in- 
puts than  the  LMS  algorithm  since  an  initial  weight  vector,  W(l),  and  a y 
value  did  not  have  to  be  chosen.  The  expected  value  algorithm  never  produced 
error  values  greater  than  100,  but  the  LMS  algorithm  produced  many  values  that  ; 
were  too  large  to  be  graphed  on  the  desired  scale. 

The  complexity  of  the  two  algorithms  is  approximately  the  same.  In 
order  to  calculate  the  weight  vector,  each  algorithm  requires  two  multipli- 
cations. The  expected  value  algorithm  also  required  an  additional  division. 
However,  in  terms  of  performance,  the  extra  division  time  is  well  worth  the 

small  trade  off  in  speed.  Each  algorithm  computes  the  same  terms,  x^x^  and 
2 

x2  , but  the  expected  value  algorithm  updates  two  estimates  at  each  Iteration 
before  dividing  the  two  to  obtain  an  estimate  of  the  weight,  W.  The  LMS  al- 
gorithm, however,  only  computes  one  estimate  of  the  weight.  Although  the  one 
estimate  takes  less  time,  the  two  estimate  approach  produces  better  perform- 
ance. 

The  problem  In  choosing  a y value  proved  to  be  great.  In  general,  as 
the  u value  decreased,  the  rate  of  convergence  also  decreased,  but  the  minimum 
values  of  the  curves  approached  more  closely  the  theoretical  minimum  error  due 
to  fewer  fluctuations  in  the  weight  vector.  Since  only  the  first  500  itera- 
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tions  were  used,  the  range  of  u was  even  smaller  than  equation  (3.3)  would 
indicate.  This  small  range  proved  to  be  an  important  drawback  to  the  IMS  al- 
gorithm. The  large  error  values  produced  by  some  values  of  y created  storage 
and  graphing  problems  that  were  not  present  In  the  expected  value  algorithm. 
For  signal -to-noise  ratios  of  1.0  and  .10,  the  reason  that  the  larger  values 
of  u generated  curves  off  the  scale  was  due  to  the  fact  that  the  y value  was 
not  in  the  range  for  convergence.  However,  for  small  values  of  li,  the  large 
errors  were  probably  due  to  the  slow  rate  of  convergence  of  the  weights. 

Although  the  choice  of  the  Initial  weight  vector  was  not  Important 
in  this  study  because  of  such  a large  number  of  iterations,  It  could  have  an 
effect  when  considering  a smaller  number.  In  any  case,  the  effect  of  the 
choice  is  not  as  Important  as  that  of  the  u value  because  it  tends  to  wash 
out  after  several  Iterations.  This  is  another  value,  however,  which  must  be 
chosen  for  the  LMS  algorithm  and  which  is  not  necessary  for  the  expected  value 
algorithm. 

The  comparison  between  the  two  algorithms  is  dramatic  when  the  slgnal- 
to-noise  ratio  Is  high,  but  becomes  even  more  so  as  the  signal -to-noise  ratio 
decreases.  When  the  signal -to-noise  rations  10.0,  the  results  are  closer  if 
a correct  y value  Is  chosen  than  when  the  slgnal-to-nolse  ratio  is  .10  and 
only  the  expected  value  algorithm  can  be  said  to  converge  within  five  hundred 
Iterations.  The  performance  for  all  signal -to-noise  ratios  and  correlation 
values  tested  reveal  that  the  expected  value  algorithm  works  better  for  a two 
Input  noise  canceling  system  with  a dc  signal  and  Gaussian  noise  than  does 
the  LMS  algorithm. 


Table  3.1.  Representative  Values  of  the  Curves  Obtained: 
Two  Input  Noise  Canceling  System 


Each  page  contains  four  graphs  labeled  a,  b,  c,  and  d. 
characteristic  values  given  by  several  parameters: 
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- the  theoretical  minimum  mean  square  error 
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Table  3.1  cont'd. 
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Figure  3.3.  Computer  Simulation  of  the  Two  Input  System. 
S/N  » 10.0  and  p12  * -75. 

(a)  Expected  value  algorithm. 

(b)  LMS  algorithm,  y * .01. 

(c)  LMS  algorithm,  y * .001. 

(d)  LMS  algorithm,  y * .00001. 
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Figure  3. 


. Computer  Simulation  of  the  Two  Input  System. 
S/N  • 10.0  and  * *50. 

(a)  Expected  value  algorithm. 

(b)  LMS  algorithm,  u * .01. 

(c)  LMS  algorithm,  w * .001. 

(d)  LMS  algorithm,  u ■ .0001. 
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Figure  3.8.  Computer  Simulation  of  the  Two  Input  System. 
S/N  * .10  and  * .99. 

(a)  Expected  value  algorithm. 

(b)  LMS  algorithm, 'ii  * .01. 

(c)  IMS  algorithm,  u * .001. 

(d)  LMS  algorithm,  u ■ .0001. 
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Figure  3.9.  Computer  Simulation  of  the  Two  Input  System. 
S/N  * .10  and  p * .98. 

(a)  Expected  value  algorithm. 

(b)  IMS  algorithm,  y * .01. 

(c)  IMS  algorithm,  y a .001. 

(d)  LMS  algorithm,  y * .0001. 
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Figure  3.10.  Computer  Simulation  of  the  Two  Input  System. 
S/N  » .10  and  * *97. 

(a)  Expected  value  algorithm. 

(b)  IMS  algorithm,  y ■ .01. 

(c)  LMS  algorithm,  u a .001. 

(d)  LMS  algorithm,  y * .0001. 


CHAPTER  IV 

EXTENSION  TO  THREE  INPUT  NOISE  CANCELING  SYSTEM: 

DERIVATION  OF  THE  ALGORITHMS 

The  algorithms  derived  and  simulated  in  Chapters  II  and  III  are  ex- 
tended in  this  chapter  to  a three  input  noise  canceling  system.  The  goal  of 
the  algorithms  is  to  minimize  the  mean  square  error.  Again  the  performance 
of  the  algorithms  is  determined  from  a computer  simulation. 

A three  input  noise  canceling  system  has  one  input  containing  a mix- 
ture of  signal  and  noise  and  two  inputs  of  noise  only.  The  weight  vector 
contains  two  weights.  The  three  input  system  is  shown  in  Figure  4.1. 


Figure  4.1.  Three  Input  Noise  Canceling  System. 

(45) 
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The  signal,  s,  is  uncorrelated  with  the  noises,  n^ , n 2,  and  n3,  but 
the  noises  are  correlated  with  each  other.  The  estimated  signal,  s,  is 
given  by 

s » x-,  - ^2*2  " W3X3  (4.1) 

and  the  error  is  defined  by  (2.2)  as 

c ■ s - s (2.2) 

The  weights,  VL  and  W^,  are  determined  such  that  the  signal  estimate, 
s,  is  a linear  least  mean  square  estimate. 


1 . LMS  Alqori thm 

The  LMS  algorithm  is  easily  extended  to  the  three  input  noise  can- 
celing system  of  Figure  4.1  because  several  scalar  quantities  of  the  two 
input  system  are  simply  replaced  with  vectors.  The  weight,  W,  of  the  two 
input  system  is  replaced  with  the  vector  W where 


The  input  Xj  is  replaced  with  a vector  of  noise  only  inputs: 


(4.2) 


1 * ft]  (4'3’ 

The  next  weight  vector,  , is  calculated  utilizing  (2.13)  and  (2.59)  with 
(4.2)  substituted  for  W and  (4.3)  substituted  for  x Thus, 

-j+1  * -J  + 2uSj~j  (4,4) 


or  in  matrix  form: 


+ 2pSj 


(4.5) 
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where 


or 


sj  ■ k 


Sj  * Xi  - W 2 x 2 ~ W«  x« 

J 'j  Zj  3j  3j 


(4.6) 

(4.7) 


In  terms  of  the  inputs,  the  weights  of  the  three  input  system  are 
given  in  Figure  4.2. 


The  IMS  algorithm  can  be  easily  extended  to  any  number  of  inputs 
necessary  to  solve  a noise  canceling  problem.  The  expected  value  algorithm 
Involves  a matrix  inversion,  but  Is  an  estimation  of  the  theoretical  optimum 
solution.  The  tradeoff  between  the  two  algorithms  is  in  terms  of  the  com- 
plexity and  performance. 
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2.  Expected  Value  Algorithm 

The  expected  value  algorithm  attempts  to  minimize  the  mean  square 
error  by  estimating  the  optimum  weight  vector  for  Figure  4.1.  The  optimum 
weight  vector  can  be  determined  by  solving  for  the  weights  that  make  the 
partial  derivative  of  the  mean  square  error  with  respect  to  the  weights 
equal  to  zero.  This  is  equivalent  to  solving  the  Wiener-Hopf  equation  for 
a three  input  system. 

Equation  (4.1)  can  be  written  in  matrix  form: 

s - x1  - WTX  (4.8) 

where  W and  X,  were  defined  in  (4.2)  and  (4.3).  Substituting  (4.2)  and  (4.3) 
into  the  definition  of  the  error  given  by  (2.2)  yields 

e * x1  - WT X - s (4.9) 

Squaring  and  taking  the  expected  value  of  both  sides  of  equation  (4.9)  pro- 
duces the  mean  square  error: 

ECe2]  * EC(x1  - s)2]  - 2WTE[(x1  - s)X]  + WTE[XXT]W  (4.10) 

In  order  to  find  the  value  of  W that  minimizes  (4.10),  it  is  necessary  to 
take  the  partial  derivative  with  respect  to  W.  Thus 

ECe2]  » -2E[(Xl  - S)XJ  + 2E[XXT]W  (4.11) 

Expanding  (4.11)  and  recalling  the  hypotheses  that  E[sx2]  * 0 and  ECsx^]  * 0 
gives 

EC*2]  - -2E[XlXj  + 2E[XXT]W  (4.12) 

Define  the  two  expected  value  terms  in  (4.12)  as  follows: 

E[XXT] 
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E[x2x2]  E[x2x3] 


— ECx-x-]  E[x,x,] 


* Etxl^ 


E[x^ x2 J 

EOVj] 


Thus,  substituting  (4.13)  and  (4.14)  into  (4.12)  gives 


• .13) 


(4.14) 


3W  E[e  ] a + 2lxx- 


(4.15) 


In  order  to  find  the  minimum  error,  the  partial  derivative  in  (4.15)  is  set 


equal  to  zero.  Solving  for  W gives 


(4.16) 


Equation  (4.16)  can  also  be  written  in  matrix  form  utilizing  (4.13)  and 


(4.14): 


E[x2x2]  E[x2x3]  W2 

E[x3x2]  E[x3x3]_  Wj 


ECx1x23 

E[x1x3] 


Taking  the  matrix  inverse,  (4.17)  can  be  solved  for  W2  and  W3: 


E[x^x2]E[x3x3]  " E[x2x3]E[x^x3] 


[x2x2JE[x3x3. 


lx2x3JE[x2x3. 


(4.17) 


(4.18) 


ECxgXgjECx^^  - E[x2x3]E[x1x2] 


.x2x2JE[x3x3J 


Lx2x3JE[x2x3. 


(4.19) 


Thus  the  values  of  W2  and  W3  as  determined  in  (4.18)  and  (4.19)  are 
the  optimum  values  necessary  to  minimize  the  mean  square  error.  They  are 
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in  fact  the  solutions  to  the  Wlener-Hopf  equation,  (1.2),  for  a three  input 
system  where  R * and  £ » ^ The  expected  value  algorithm  utilizes 
these  optimum  weights  in  order  to  estimate  the  signal.  The  expected  values 
are  mathematical  terms  and  are  estimated  (see  (2.61))  for  the  computer  simu- 
lation of  the  three  input  system.  Figure  4.3  shows  the  three  Input  noise 
canceling  system  with  the  weights  determined  from  the  theoretical  derivation 
of  the  expected  value  algorithm. 
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Since  the  expected  value  and  LMS  algorithms  both  converge  to  the  same  error 
(see  Appendix  A),  It  Is  useful  to  determine  the  theoretical  minimum  mean 
square  error. 


3.  Theoretical  Minimum  Mean  Square  Error 

In  order  to  effectively  compare  the  two  algorithms,  it  is  desirable 
to  determine  the  theoretical  minimum  mean  square  error.  The  optimum  weight 
vector  for  a linear  least  mean  square  estimate  was  determined  by  the  expected 
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value  algorithm.  Thus  using  (4.18)  and  (4.19)  for  the  weight  vector,  (4.1) 
for  the  estimated  signal,  and  (2.40)  for  the  minimum  error,  the  theoretical 
minimum  error  can  be  calculated: 

'min  ' E«s'  ' S>21  <2-40> 
and 


eMMS  * E^nl  ~ W2X2  " W3X3^^  (4.20) 

Equation  (4.20)  can  be  simplified  using  Figure  4.1  and  replacing  , x^,  and 
x^  with  their  equivalent  values  in  terms  of  s,  n-j , n2,  and  n3.  The  definition 
of  correlation  given  by  (2.34)  Is  also  used: 

eMMS  * ECnl2]  • 2w2ECnln2]  ‘ 2W3EEV3]  + 2W2W3E[n2n3] 

+ W22E[n22]  + W32E[n32] 

and  thus 

eMMS  * E^nl  ^ " 2^2p1  Etn,]2]Etn2^  “ 2^3p1  3 VE^nl ^lE^n3^] 

+ 2W2W3P23VECn22]ECn32^  + W22ECn22]  + (4-21) 


The  optimum  weight  vector  given  by  (4.18)  and  (4.19)  can  be  rewritten 
in  terms  of  s,  n^ , n2«  and  n3-  The  definition  of  correlation  in  (2.34)  also 
helps  to  simplify  the  vector.  Thus, 


ft"  1 : 

EC",2] 


(4.22) 


where 


12 


' p23p13 

T 


1 - p 


23 


(4.23) 
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and 


W * c-i 

3 2 


f«"i  ] 
Et»,Z] 


where 


p - p p 
. _ 13  12  23 

C2  ' { 2 

1 ‘ p23 


(4.24) 


(4.25) 


Substituting  (4.22)  and  (4.24)  into  (4.21),  the  theoretical  minimum 
error  can  be  simplified  further: 


’MMS 


* 2^-j  * * ^1^2^23  ^ ^1  ^ ^2  ^ (4*26) 


where  c-|  and  are  defined  by  (4.23)  and  (4.25). 

Since  the  weight  vectors  of  both  algorithms  converge  to  the  optimum 
weight  vector  (see  Appendix  A),  the  theoretical  minimum  mean  square  error  of 
both  algorithms  is  given  by  (4.26). 


4.  Summary  of  the  Algorithms 

The  algorithms  are  sumnarized  here  in  the  forms  that  are  used  for  the 
computer  simulation.  The  inputs,  weights,  and  estimates  are  presented  In 
terms  of  the  jth  iteration. 

A.  LMS  Algorithm 

The  LMS  algorithm  utilized  the  computed  value  of  Sy^  for  the  com- 
puter simulation.  The  values  of  the  weights  at  each  iteration  are  stored 
so  that  they  can  be  used  to  compute  the  weights  for -the  next  iteration. 
Therefore  the  IMS  algorithm  is  simulated  for  the  jth  iteration  as  shown  In 
Figure  4.4. 
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Figure  4.4.  LMS  Algorithm  as  Computer  Simulated. 


B.  Expected  Value  Algorithm 

The  expected  value  algorithm  is  derived  in  terms  of  the  expected 
values  of  the  inputs  which  will  minimize  the  mean  square  error.  Since  the 
expected  values  are  mathematical  quantities  which  are  not  known  exactly,  the 
expected  value  algorithms  must  estimate  them  and  is  therefore  only  an  approxi- 
mation to  the  optimum  solution. 

The  expected  values  were  estimated  using  the  iterative  procedure  des- 
cribed in  Chapter  II.  The  expected  value  algorithm  for  the  jth  iteration  of 
the  computer  simulation  is  implemented  as  seen  in  Figure  4.5. 


C.  Comparison  of  the  Algorithms 

The  optimum  weight  vector  will  produce  the  best  results  for  the  three 
Input  system.  However,  since  It  can't  be  calculated  directly  the  expected 
value  algorithm  approximates  It  by  estimating  the  expected  values.  As  in 
the  two  Input  case,  the  LMS  algorithm  attempts  to  estimate  the  optimum  weight 
vector  as  one  term  while  the  expected  value  algorithm  estimates  each  expected 
value  before  combining  them  to  form  the  weight  vector.  Although  this  method 
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Figure  4.5.  Expected  Value  Algorithm  as  Computer  Simulated. 

is  more  complex  and  more  time  consuming  due  to  the  increased  number  of  mul- 
tiplications and  divisions,  the  improved  performance  justifies  It.  The  U1S 
algorithm  could  also  be  written  in  terms  of  a feedback  system  as  in  the  two 
Input  case,  but  this  is  not  necessary.  The  expected  value  algorithm  can't 
be  a feedback  system  because  it  needs  the  jth  inputs  to  compute  the  jth  es- 
timate. The  LMS  algorithm  instead  computes  the  jth  weights  with  the  data 
through  the  (j-l)st  inputs.  This  method  is  not  as  effective  as  that  of  the 
expected  value  algorithm. 


CHAPTER  V 

COMPUTER  SIMULATION  OF  THE  THREE  INPUT 
NOISE  CANCELING  SYSTEM 

The  three  input  algorithms  developed  in  Chapter  IV  were  simulated  on 
a PDP-11/45  computer  at  Duke  University.  Figures  4.4  and  4.5  give  the  al- 
gorithms in  terms  of  their  computer  implementation. 

1 . Problem  Description 

Parameter  values  had  to  be  chosen  before  starting  the  simulation. 

The  signal  was  a dc  signal  with  value  10.0  as  in  the  two  input  case.  The 
values  of  the  initial  weight  vector  for  the  LMS  algorithm  were  again  constant. 
The  first  weight,  Wg,  was  set  equal  to  5.0  as  in  the  two  input  case,  and  the 
second,  W^,  was  arbitrarily  set  to  3.0.  These  values  were  not  important  due 
to  the  large  number  of  iterations. 

The  choice  of  the  parameter  » for  the  LMS  algorithm  was  again  a 
problem  as  It  was  for  the  two  input  system.  In  order  for  the  algorithm  to 
converge,  the  value  of  u must  be  in  the  range  (see  Appendix  A)  defined  as: 

|1  - 2uE[XTX]|  < 1 (5.1) 

where  £ is  defined  by  (4.3).  The  three  u values  of  .01,  .001,  and  .0001, 
corresponding  to  graphs  b,  c,  and  d,  respectively  on  each  page  of  figures, 
are  the  same  as  those  used  In  the  two  Input  system  simulation.  Also  for 
convenience  in  the  simulation. 
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(5.2) 


St",2]  ■ E[n22]  * £[n32] 

When  the  signal-to-noise  ratio  is  10.0,  all  three  values  of  u are  within  the 
range  for  convergence  of  the  algorithm.  When  the  signal-to-noise  ratio  is 
lowered  to  1.0,  the  u value  of  .01  is  outside  the  range  of  convergence. 

Both  u values  of  .01  and  .001  are  outside  the  range  of  convergence  when  the 
signal-to-noise  ratio  is  .1.  For  a given  signal-to-noise  ratio,  a decrease 
in  the  parameter  u results  In  a decrease  in  the  rate  of  convergence  since 
the  quantity  |1-2yE[X/Y]  [ approaches  one. 

The  noise  was  generated  using  the  method  described  in  Appendix  B. 
Independent  zero  mean  Gaussian  random  numbers  were  generated  and  then  trans- 
formed using  this  method  to  correlated  Gaussian  random  numbers.  The  correla- 
tions, p 12,  p^,  an<*  Pg 3’  were  Parameter  values  of  the  program.  The  input 
matrix  E[X.xJ]  has  to  be  a positive  semi -definite  matrix.  Thus  the  trans- 

J J 

formation  matrix,  made  up  of  correlation  values  and  given  by  (5.3) 


(5.3) 

has  to  positive  semi-definite.  This  restricts  the  choices  of  the  correlation 
values.  Once  two  correlation  values  have  been  chosen,  the  third  must  be  se- 
lected in  such  a manner  that  the  matrix  is  positive  semi -definite. 

The  mean  square  error  was  approximated  using  (3.2)  as  In  the  two  In- 
put case.  Iterative  estimates  were  used  for  the  expected  values  as  shown  in 
Figure  4.5. 

The  theoretical  minimum  mean  square  error  was  defined  in  (4.26).  Its 

2 

value  depends  only  upon  the  values  of  E[n^  ],  p^»  0^3.  and  033-  when  the 
three  correlations  all  have  values  of  zero,  the  worst  case  is  realized  and 
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no  noise  can  be  effectively  canceled.  As  the  correlations  increase  in  value, 


noise  canceling  occurs. 

The  method  of  comparison  of  the  two  algorithms  was  the  same  as  for 
the  two  input  case.  Graphs  of  the  mean  square  error  as  a function  of  the 
number  of  iterations  were  used.  Twenty-five,  five  hundred  iteration  ensem- 
bles were  averaged  to  obtain  the  curve  that  was  plotted.  As  in  the  two  input 
system,  the  three  input  system  simulation  curves  are  compared  for  differences 
in  rates  of  convergence  and  minimum  error  after  500  iterations. 

2.  Performance  Comparison 

Following  are  the  graphs  of  the  mean  square  error  as  a function  of 
the  number  of  iterations.  Each  page  contains  graphs  with  the  same  signal- 
to-noise  ratio  and  correlation  values.  One  graph  using  the  expected  value 
algorithm  and  graphs  of  three  values  of  u for  the  LMS  algorithm  are  on  each 
page.  Again  all  of  the  graphs  are  plotted  on  the  same  scale.  The  straight 
line  effect  (see  Appendix  C)  is  still  observed  when  the  value  of  the  mean 
square  error  is  too  large  for  the  boundaries  of  the  scale.  Table  5.1  gives 
representative  values  for  the  curves. 

Figures  5.1  through  5.3  set  the  signal-to-noise  ratio  to  10.0. 

I 

Each  figure  has  a different  set  of  correlation  values.  The  expected  value 
algorithm  exhibits  the  fastest  rate  of  convergence  and  the  smallest  error 
after  five  hundred  Iterations.  However,  graphs  5.1b  and  5.2b,  which  are 
the  u values  of  .01  for  the  LMS  algorithm,  are  very  similar  to  those  of  the 
expected  value  algorithm.  Thus  if  a "correct"  u value  Is  chosen,  the  re- 
sults can  be  close  to  the  optimum  results  as  estimated  by  the  expected  value 
algorithm.  As  the  u value  decreases,  a slower  rate  of  convergence  is  evi- 
denced and  the  minimum  error  after  five  hundred  iterations  is  further  from 
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the  theoretical  minimum.  If  two  of  the  correlation  values  are  small  (Fig- 
ure 5.3),  the  difference  between  the  expected  value  and  LMS  algorithm  be- 
comes much  greater. 

Figures  5.4  through  5.6  lower  the  signal-to-noise  ratio  to  1.0.  The 
same  correlation  values  are  used  again.  The  smaller  signal-to-noise  ratio 
produces  dramatic  results.  The  expected  value  algorithm  generates  the  best 
curves  in  terms  of  the  rate  of  convergence  and  the  minimum  error.  The  v 
value  of  .001  for  the  LMS  algorithm  produces  the  only  curve  that  can  truly 
be  compared  with  the  expected  value  algorithm.  Nevertheless,  its  rate  of 
convergence  is  much  slower  and  its  minimum  error  is  much  greater.  The  u val- 
ue of  .01  generally  produces  values  that  are  too  large  for  the  desired 
scale  or  curves  that  appear  to  converge  to  a large  error.  For  this  signal- 
to-noise  ratio,  the  u value  of  .01  is  outside  the  range  necessary  for  con- 
vergence. The  u value  of  .0001  has  such  a slow  rate  of  convergence  that  the 
straight  line  effect  is  again  observed.  As  the  correlation  decreases,  this 
case  of  the  LMS  algorithm  has  an  even  slower  rate  of  convergence  and  even- 
tually the  value  of  the  mean  square  error  for  five  hundred  iterations  is  too 
large  for  the  desired  scale. 

Figures  5.7  through  5.9  lower  the  signal-to-noise  ratio  to  .10.  Only 
the  expected  value  algorithm  generates  curves  with  values  small  enough  for 
the  desired  scale.  The  u values  of  .01  and  .001  are  outside  the  range  for 
convergence  of  the  algorithm.  When  u equals  .0001,  the  rate  of  convergence 
is  too  slow  for  convergence  within  five  hundred  iterations  and  the  minimum 
values  of  the  curves  are  still  too  large  for  the  scale.  Although  the  u value 
of  .00001  Is  not  shown,  the  minimum  value  was  found  to  be  larger  (for  only 
five  hundred  Iterations)  than  the  case  where  y equaled  .0001.  Again  there 
is  a storage  and  graphing  problem  created  due  to  the  large  numbers  generated 
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by  the  LMS  algorithm  that  is  not  present  with  the  expected  value  algorithm. 

3.  Summary  of  the  Three  Input  System  Results 

The  expected  value  algorithm  proved  to  be  the  superior  algorithm 
for  the  three  Input  case  as  it  was  for  the  two  input  system.  The  rate  of 
convergence  is  consistently  faster  and  the  minimum  error  is  smaller  than  for 
any  of  the  u values  of  the  LMS  algorithm  that  were  tested.  For  a signal  - 

! i 

to-noise  ratio  of  10.0,  the  LMS  algorithm  produces  curves  similar  to  those 
of  the  expected  value  algorithm.  However,  this  occurs  for  only  one  particu- 
lar u value  since  other  y values  do  not  produce  results  that  are  as  similar. 

As  the  signal-to-noise  ratio  decreases,  the  expected  value  algorithm  clearly 
produces  the  best  results.  Several  y values  also  generate  curves  which 
have  values  too  large  to  be  plotted. 

In  several  cases,  the  expected  value  algorithm  generated  curves 
whose  minimum  values  were  less  than  the  theoretical  minimum  values.  This  Is 
probably  a peculiarity  of  the  computer  simulation  and  the  computer  used. 

Since  all  of  the  noise  was  generated  with  the  same  variance,  this  probably 
contributed  to  the  small  minimum  errors.  This  would  be  especially  noticeable 
with  small  signal-to-noise  ratios  because  of  a higher  possibility  of  arith- 
metic roundoff  errors. 

The  expected  value  algorithm  is  more  complex  in  terms  of  the  number 
of  multiplications  and  divisions  required  to  implement  it.  Again,  the 
several  estimate  approach  of  the  weight  vector  used  by  the  expected  value 
algorithm  enables  a better  estimate  in  fewer  iterations.  For  a three  Input 
noise  canceling  system,  a two  dimensional  matrix  Inverse  Is  required.  Al- 
though an  inverse  of  this  size  matrix  Is  more  complex  than  the  LMS  algorithm, 
the  tradeoff  in  terms  of  performance  appears  to  justify  this  added  complexity. 
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The  effect  of  changes  In  correlation  Is  hard  to  assertain  because  of 
the  large  number  of  values  per  curve.  Generally  as  the  correlations  were 
lowered,  the  difference  between  performance  of  the  two  algorithms  became 
greater.  With  three  inputs,  the  signal -to-noise  ratio  appeared  to  have  more 
of  an  effect  on  the  results. 

The  problem  with  the  LMS  algorithm  is  the  same  as  in  the  two  Input 
system.  Values  for  the  Initial  weight  vectors  and  for  u must  be  chosen. 

It  Is  almost  impossible  to  determine  beforehand  the  best  values  for  these 
parameters  (in  terms  of  minimum  mean  square  error  and  convergence  rate).  Un- 
fortunately, unless  the  u value  is  chosen  correctly,  the  LMS  algorithm  may 
not  converge  or  may  have  a slow  rate  of  convergence  and  large  errors  even 
after  five  hundred  iterations. 


! 
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Table  5.1.  Representative  Values  of  the  Curves  Obtained: 
Three  Input  Noise  Canceling  System 


Each  page  contains  four  graphs  labeled  a,  b,  c,  and  d.  These  have 
characteristic  values  given  by  several  parameters: 

S 

N - the  signal-to-noise  ratio 

p-|2  " the  correlation  between  the  noises  and 

p-|3  - the  correlation  between  the  noises  n-j  and  n^ 

?23  - the  correlation  between  the  noises  n^  and  n^ 

TMIN  - the  theoretical  minimum  mean  square  error 

Algorithm  - the  algorithm  employed-either  expected  value  or  IMS 

u - the  value  of  the  parameter  » for  the  LMS  algorithm 

MAX  - the  maximum  value  of  the  curve 

MIN  - the  minimum  value  of  the  curve 
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Table  5.1  cont'd 
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Figure  5.6.  Computer  Simulation  of  the  Three  Input  System. 

S/N  * 1.0,  p.j2  * -25,  p^2  * *99,  and  P22  r *25. 

(a)  Expected  value  algorithm. 

(b)  IMS  algorithm,  u * .01. 

(c)  IMS  algorithm,  p « .001. 

(d)  LMS  algorithm,  p = .0001. 
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Figure  5.7.  Computer  Simulation  of  the  Three  Input  System. 

$/N  ■ .10,  p.|2  * -90,  p13  * .90,  and  p23  * .80. 

(a)  Expected  value  algorithm. 

(b)  LMS  algorithm,  u * .01. 

(c)  LWS  algorithm,  u ■ .001. 

(d)  LMS  algorithm,  u * .0001. 
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Figure  5.8.  Computer  Simulation  of  the  Three  Input  Systan. 

S/N  * .10,  * -75,  ^ * -20,  and  P23  * *75. 

(a)  Expected  value  algorithm. 

(b)  LMS  algorithm,  u * .01. 

(c)  LMS  algorithm,  u * .001. 

(d)  LMS  algorithm,  u * .0001. 
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Figure  5.9.  Computer  Simulation  of  the  Three  Input  System. 

S/N  * .10,  p12  * *25,  * .99,  and  p23  « .25. 

(a)  Expected  value  algorithm. 

(b)  IMS  algorithm,  y * .01. 

(c)  LMS  algorithm,  y * .001. 

(d)  LMS  algorithm,  y a .0001. 


CHAPTER  VI 

SUMMARY  AMD  FUTURE  WORK 

1 . Summary 

This  research  compares  two  algorithms  used  in  noise  canceling  systems. 
The  first  part  studies  a two  input  noise  canceling  system.  The  expected  val- 
ue and  LMS  (least  mean  square)  algorithms  are  compared  in  terms  of  their  mean 
square  error  performance.  The  second  section  extends  both  algorithms  to  a 
three  input  noise  canceling  system.  These  results  are  then  compared. 

Chapter  II  derives  the  algorithms  for  the  two  input  noise  canceling 
system.  The  optimum  weight  vector  is  given  by  the  Wiener-Hopf  equation  as 
explained  in  Chapter  I.  The  expected  value  algorithm  as  implemented  for  the 
computer  simulation  is  an  estimate  to  the  optimum  algorithm  in  that  the  math- 
ematical expected  values  are  replaced  with  estimates.  The  LMS  algorithm 
attempts  to  lessen  the  complexity  involved  in  the  matrix  inversion  by  using 
gradient  techniques  and  the  method  of  steepest  descent. 

Chapter  III  reveals  the  results  of  the  two  input  system  computer 
simulation.  The  expected  value  algorithm  proved  superior  in  that  it  had  a 
faster  rate  of  convergence  and  lower  minimum  errors  after  five  hundred  iter- 
ations. The  difference  was  especially  apparent  at  lower  slgnal-to-nolse 
ratios  and  lower  values  of  noise  correlation.  The  LMS  algorithm  also  pre- 
sented problems  in  that  values  for  the  Initial  weight  vector  and  for  u had  to 
be  chosen.  Although  the  choice  of  the  initial  weight  vector  was  not  as 
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important  when  considering  five  hundred  iterations,  the  value  of  u was  ex- 
tremely Important.  In  order  for  the  algorithm  to  converge,  the  value  of  the 
term  had  to  be  within  a specified  range  which  proved  small  when  only  five 
hundred  iterations  were  considered. 

The  algorithms  are  extended  to  a three  input  noise  canceling  system  in 
Chapter  IY.  Again  the  expected  value  algorithm  is  an  estimate  to  the  optimum 
for  a minimum  linear  least  mean  square  estimate.  The  IMS  algorithm  is  much 
simpler  because  it  does  not  involve  a matrix  inversion. 

Chapter  V reveals  the  results  of  the  computer  simulation  of  the  three 
input  system.  The  results  are  similar  to  those  for  the  two  input  system. 

The  expected  value  algorithm  has  a faster  rate  of  convergence  and  smaller 
minimum  errors.  However,  for  several  cases  of  correlation  values,  the  LMS 
algorithm  produced  results  that  were  extremely  similar  to  those  produced  by 
the  expected  value  algorithm.  Nevertheless,  this  only  happened  for  one  value 
of  u while  the  other  two  u values  tested  produced  results  that  were  much  worse. 
The  additional  complexity  of  the  expected  value  algorithm  due  to  the  matrix 
Inverse  is  still  justified.  Unless  the  "correct"  u value  is  chosen,  the  re- 
suits  will  be  much  'worse  than  the  optimum. 

In  general,  the  expected  value  algorithm  performs  better  than  the  LMS 
algorithm.  For  a noise  canceling  system  with  only  cwo  or  three  inputs,  an 
estimate  to  the  optimum  method  is  better.  With  a small  number  of  inputs,  the 
added  complexity  involved  in  Inverting  the  input  matrix  Is  justified  in  terms 
of  performance.  This  added  complexity  in  forming  the  weight  vector  is  prob- 
ably the  reason  for  the  better  performance.  Although  both  algorithms  compute 
the  same  terms  (in  terms  of  the  data),  the  expected  value  algorithm  updates 
several  estimates  before  combining  these  to  form  one  estimate  for  the  weight 
vector.  The  LMS  algorithm,  however,  only  forms  one  estimate.  The  several 
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estimate  approach  of  the  expected  value  algorithm  produces  a better  estimate 
of  the  weight  vector  in  a fewer  number  of  iterations. 

The  additional  problems  encountered  when  implementing  the  LMS  algor- 
ithm (that  is,  choosing  values  for  u and  the  initial  weight  vector)  also 
discourage  its  use.  The  expected  value  algorithm  does  not  have  these  prob- 
lems and  its  performance  is  better  in  terms  of  rate  of  convergence  and  min- 
imum errors. 

2.  Future  Research 

This  research  should  encourage  efforts  to  study  the  tradeoffs  between 
complexity  and  performance  that  exist  in  many  algorithms.  Often  algorithms 
are  developed  for  the  more  intricate  cases  of  a problem.  Thus  while  the 
algorithm  may  be  the  only  practical  way  to  solve  the  more  complex  problem, 
this  is  not  always  true  when  the  problem  is  in  a simpler  form. 

The  two  input  noise  canceling  system  was  studied  thoroughly  and  the 
expected  value  algorithm  produced  better  results.  The  three  input  system 
needs  to  be  studied  more  carefully.  The  relationship  between  changes  in 
correlation  values  should  be  established.  With  variations  for  three  correla- 
tion values,  only  a few  were  selected  to  study  here.  The  Importance  of  one 
correlation  value  with  respect  to  the  others  should  be  ascertained.  Changes 
In  the  variance  between  the  noise  Inputs  should  also  be  studied. 

In  general,  the  importance  of  the  initial  weight  vector  of  the  LMS 
algorithm  should  be  studied.  For  this  research,  the  value  was  arbitrarily 
chosen.  If  less  than  five  hundred  Iterations  are  essential,  the  Initial 
weight  value  will  be  more  Important. 

This  study  should  also  be  extended  to  systems  with  more  Inputs  to 
determine  how  large  a system  can  be  handled  practically  with  the  expected 
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value  algorithm.  Many  programs  today  are  available  to  invert  a matrix  quick- 
ly. This  would  enable  larger  systems  to  employ  this  estimate  to  the  optimum 
weight  vector. 

This  research  is  only  a beginning,  but  It  should  encourage  further 
research  in  this  area.  The  major  contribution  should  be  the  questioning  of 
algorithms  to  solve  simpler  problems  when  the  optimum  method  would  produce 
much  better  results. 


APPENDIX  A 


CONVERGENCE  OF  THE  WEIGHT  VECTOR  FOR 
THE  LMS  ALGORITHM 


In  order  to  make  the  comparison  between  the  expected  value  and  LMS 
algorithms  valid,  it  is  necessary  that  their  weight  vectors  converge  to  the 
same  value.  The  weight  vector  for  the  expected  value  algorithm  is  shown  to 
be  an  estimate  to  the  optimum  weight  vector  in  Section  3 of  Chapter  II. 
Therefore,  it  is  only  necessary  to  prove  that  the  weight  vector  for  the  LMS 
algorithm  converges  to  the  optimum  weight  vector. 

Looking  first  at  the  two  input  noise  canceling  system,  the  weight 
vector  for  the  LMS  algorithm  is  defined  to  be 

W, * W,  + 2ys_.x„  (A.l) 


nk  * mc  2ys  .x_ 
LMSj+l  LMSj  J 2j 


The  optimum  and  theoretical  expected  value  algorithms  have  a weight  vector 
given  by 


E[x,  x«  ] 

W . J J 

<»'}  evj  ieCy 


(A. 2) 


To  prove  that  the  LMS  weight  vector  converges  to  the  optimum  weight 


vector,  it  is  essential  to  show  that 
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For  the  two  input  noise  canceling  system, 

' “jX2j 

Substituting  (A.l)  and  (A.4)  into  (A. 3)  produces 


EtWLMSj+]^  * EC(1  " 2uxZ.  ^MLMSj  + 2uXTjX2j‘! 


(A.4) 


(A. 5) 


Reducing  the  index  by  one,  can  be  written  in  terms  of  Wy^ . • Thus 

(A. 5)  becomes 


E[Wlms  ] - E[(l-2ux,  2)(1-2ux2  2)W.M5  + (1-2ux  2)(2ux,  x ) 

LMSj+l  2j  2j-l  ^j-l  2j  2j-1 


+ 2ux,  x-  ] 

'j  Zj 


(A. 6) 


Equation  (A. 6)  can  be  rewritten  in  terms  of  Both  inputs  x-j  and  x2  are 

uncorrelated  over  time.  Therefore,  the  expected  value  can  be  taken  of  each 
product  in  (A. 6): 

E^WLMSj+1^  * C1 ”2u£Cx22]) J 1wlmsq  + 

2uE[x1x2][1  + (1-2uE[x22])  + (1-2uE[x22])2  + ... 

+ (l-2uE[x22])j]  (A. 7) 


Recognizing  that  the  second  term  in  (A. 7)  is  a geometric  series,  (A. 7)  can 


be  simplified: 


2 1+1  l-Cl-2uE[x  2]]J+1 

(l-2yE[x \ 2vE[Xlx2]  (A. 8) 

2 ^0  _1-[l-2vE[x  2]]  . 1 2 


In  order  for  (A. 8)  to  converge,  it  is  required  that 

|1  - 2uE[x22]|  < 1 


(A. 9) 


Assuming  that  (A. 9)  is  true,  taking  the  limit  as  j approaches  infinity  of 
both  sides  of  (A. 8)  gives 
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lim 
j -*■  « 


E[“usJ+1J 


E[x, 1 
E[x22] 


(A. 10) 


This  is  the  desired  result.  The  weight  vector  for  the  LMS  algorithm 
converges  to  the  optimum  weight  vector  when  (A. 9)  is  true. 

A study  of  convergence  was  made  in  more  detail  with  an  arbitrary  set 
of  parameters  for  the  two  input  noise  canceling  system.  A signal -to-noise 
ratio  of  10.0  and  a correlation  of  .95  were  used  for  the  study.  The  values 
of  u tested  were  .01,  .001  , .0001,  and  .00001.  Figure  A.l  shows  tl,e  graphs 
of  the  mean  square  error  versus  the  number  of  iterations  for  one  ensemble 
and  the  above  parameter  values.  One  ensemble  was  graphed  instead  of  the 
usual  average  of  twenty-five  ensembles  to  emphasize  the  fluctuations  in  the 
error.  Looking  at  Figure  A.l,  the  largest  u value  coincides  with  the  fastest 
rate  of  convergence,  but  the  mean  square  error  does  not  steadily  or  smoothly 
decrease  to  a minimum  value.  Instead  the  value  fluctuates  considerably. 

The  second  graph  with  u a .001  exhibits  the  smoothest  decline  towards  a min- 
imum value,  but  the  rate  of  convergence  is  slower  than  the  curve  in  (a)  and 
the  mean  square  error  still  varies  with  a small  number  of  iterations.  Never- 
theless, the  mean  square  error  in  (b)  has  reached  a smaller  value  after  five 
hundred  iterations  than  in  (a).  Graph  (c)  displays  an  even  slower  rate  of 
convergence  than  (b)  and  graph  (d)  does  not  appear  to  begin  to  converge  after 
five  hundred  iterations. 

For  the  above  parameters,  the  value  of  the  optimum  weight  is  approxi- 
mately .95.  Examining  the  value  of  the  weight  vector  at  each  iteration  re- 
veals more  about  the  rate  of  convergence  towards  the  optimum  weight  vector. 
Although  with  u ■ .01,  the  rate  of  convergence  Is  the  fastest,  there  are 
large  fluctuations  in  the  value  of  the  weight  vector.  Even  after  two  thousand 
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Iterations,  the  weight  vector  does  not  appear  to  converge,  but  instead  fluc- 
tuates in  value  between  -1  and  +3.  With  y * .001,  the  rate  of  convergence 
is  slower,  but  the  weight  vector  exhibits  fewer  variations  in  value.  The 
fluctuations  are  mainly  between  .7  and  1.7.  The  case  with  y * .0001  shows 
an  even  slower  rate  of  convergence  than  the  previous  case,  but  the  fluctua- 
tions in  the  values  of  the  weight  vector  have  almost  disappeared.  Inspecting 
the  rate  of  convergence  for  two  thousand  iterations  indicates  that  the 
weight  should  converge  after  approximately  twenty-five  hundred  iterations. 

The  last  case  studied  was  with  u * .00001.  The  weight  vector  again  showed 
almost  no  fluctuations  in  value,  but  Instead,  a steady  decline  towards  the 
optimum  value.  Convergence  based  on  this  rate  will  not  occur  until  approxi- 
mately six  thousand  iterations. 

The  results  of  the  convergence  test  can  be  extended  to  the  three  in- 
put case.  Although  the  study  is  not  as  extensive  as  for  the  two  input  sys- 
tem, it  is  shown  here  that  the  weight  vector  for  the  LMS  algorithm  converges 
to  the  optimum  weight  vector.  It  was  shown  in  Section  2 of  Chapter  IV  that 
the  weight  vector  for  the  expected  value  algorithm  is  an  estimate  to  the 
optimum  weight  vector. 

The  weight  vector  for  the  IMS  algorithm  can  be  written  in  vector 

form: 

+ 2usj-j  (AJ1) 


or 


“IMS 


j+1 


LW3. 


’“IMS, 


+ 2ys . 
0 


(A. 12) 


It  is  necessary  that 
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I. 


j1!1" . * Wopt 

The  optimum  weight  vector  is  given  by  (4.17)  and  is 

-1 


Sopt  ' %x 


or 


E[x2x2] 

E[x2x3]‘ 

-1 

ECx-,x2] 

_E[x2x3] 

E[x3x3]_ 

_E[x1x3]_ 

W 

~opt 


Taking  the  expected  value  of  both  sides  of  (A. 11)  gives 


(A.l 3) 


(A. 14) 


(A. 15) 


EtyLMSj+1]  3 E^MSj  + 2uSjV 


Substituting  for  s.  in  (A. 16)  yields 

J 


(A. 16) 


(A. 1 7) 


E^LMSj+1  ^ * E^2wXlJ.-j  + ^LMSj  ^ "2u— jT— j ^ 

Equation  (A. 17)  can  be  rewritten  in  terms  of  the  (j-l)st  iteration: 

Equation  (A.  18)  can  also  be  written  in  terms  of  Wy^  and  the  expected  value 
can  be  taken  of  each  term.  Thus 

E&us  ] - Wn(l-2uE[XTXj)j+1  + 2uE[x,Xj[l  + (1-2uE[XTX)]  + ... 
j+1 

+ (1-2vE[XTX])j]  (A. 19) 

Recognizing  that  (A. 19)  contains  a geometric  series  which  will  converge  when 

|1  - 2uE[XTX]|  < 1 (A. 20) 

(A. 19)  can  be  written  in  a simpler  form: 
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E[Wlms  ] - W0(l-2uE[XTX])j+1  + EuECx.Xj  (A.21) 

LM5j+l  0 " 1 l-[l-2uE[XTXj] 

Thus  when  (A. 20)  holds,  the  series  will  converge.  Taking  the  limit  of  both 

sides  of  (A.21)  produces 

/l".  ECWub  ] ■ (E[XTa)',E[x]a  (A. 22) 

J ' 

This  is  the  desired  result.  Thus  the  weight  vector  of  the  LMS  algor 
ithm  does  in  fact  converge  to  the  optimum  weight  vector  for  both  the  two  and 
three  input  noise  canceling  systems. 


«fisn  squire  erro 


APPENDIX  B 

NOISE  GENERATION  FOR  THE  THREE  INPUT  NOISE 
CANCELING  SYSTEM 

The  computer  simulation  requires  three  correlated  Gaussian  random 
numbers  for  the  noise  components  to  the  three  input  noise  canceling  system. 
An  unpublished  paper  by  Charles  S.  Liu,  University  of  Michigan-Dearbom, 
develops  a technique  to  transform  a n-dimensional  vector  of  independent 
Gaussian  random  numbers  to  correlated  Gaussian  random  numbers.  The  method 
will  be  summarized  here. 

Let  zj 3 [Z^ , Z2,. . . ,Zn]  be  an  n-dimensional  random  vector  whose  ele- 
ments are  independent  Gaussian  random  numbers  with  zero  mean  and  unit  var- 
iance. The  vector  1 may  be  transformed  such  that  it  contains  new  random 
variables: 

r 1 - pz  (b.i) 

where  Ms  an  n x n transformation  matrix  specified  by  the  covariance  matrix 
of  Z ’ : 

SL  ■ E[Z'Z,T] 

£«PPT  (B.2) 

The  transformation  Z'  * PZ  Is  not  unique  since  there  are  an  infinite  number 
of  matrices  which  satisfy  (B.2).  Since  the  covariance  matrix  is  real 
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symmetric,  the  eigenvectors,  (Y^}  of  £ are  orthogonal.  Thus  £ can  be  diagon- 
alized as  shown: 


CL]  »Ig » • • • » Ig » • • • » Yfll 


'1 


(B.3) 


IF  (8.3)  is  premultiplied  by  CX.^  »Xg » - - - »Xn3  and  postmulti plied  by 
CL,,^,...,^]1,  then  ' 

[Yq,Y2,...,Yfl][L1,Y2,...,Yfi]TgiCYq,Y2,...Yfi][Y^,Y2,...,Yti]T  - £ 


0, 


X1  0 


[Y_]  ‘Igt • • • 


(B.4) 


Thus  Q.  can  be  rewritten  using  (3.4)  so  that 


a - C,.le.--y 


V*T  °~ 

^ 0 

2. 

S f^2 

0 

o 'jr 

Opij yT 


(B.5) 


Recognizing  that  (8.5)  is  of  the  same  form  as  (B.2),  the  transformation 
matrix  P can  be  defined  as 


* • a,.** y 


2. 

o 


(8.6) 


Therefore  the  eigenvalues  and  eigenvectors  of  the  covariance  matrix 
determine  the  transformation  matrix  necessary  to  obtain  correlated  Gaussian 
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random  numbers  from  independent  ones.  That  transformation  is  defined  by  (B.l) 
where  P_,  the  transformation  matrix,  is  given  by  (B.6). 


APPENDIX  C 

STRAIGHT  LINE  EFFECT  OF  GRAPHS 


In  order  to  make  the  comparison  between  the  algorithms  clearer  and 
more  apparent,  it  was  desirable  to  plot  all  of  the  mean  square  error  curves 
using  the  same  graph  scale.  The  graphs  were  produced  on  a computer  screen 
of  limited  size.  When  the  values  were  too  large  to  be  within  the  range  of 
the  scale,  the  curve  was  plotted  as  a straight  line  at  the  top  of  the  screen. 
Although  the  straight  line  section  of  the  curve  tends  to  Imply  that  the 
curve  has  a constant  value  over  this  range,  this  is  not  the  case.  Figure 
C.l  shows  this  effect  by  rescaling  one  graph  so  that  the  entire  curve  lies 
within  the  boundaries  of  the  scale.  This  is  done  using  the  maximum  and  min- 
imum values  of  the  curve  as  the  maximum  and  minimum  points  of  the  scale.  If 
a graph  is  replaced  by  the  statement  that  the  minimum  error  Is  too  large  for 
the  scale,  this  means  that  for  five  hundred  iterations  the  error  Is  always 
greater  than  140. 
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APPENDIX  D 

COMPUTER  SIMULATION  PROGRAMS 

A.  Two  Input  Noise  Canceling  System 
Purpose : 

The  computer  simulation  was  the  basis  of  comparison  between  the  ex- 
pected value  and  LMS  algorithms.  The  equations  given  for  the  algorithms  in 
Figures  2.12  and  2.13  were  implemented  on  a PDP-11/45  at  Duke  University. 
The  mean  square  error  was  computed  and  plotted  as  a function  of  the  number 
of  Iterations. 

Specifications: 

The  specifications  of  the  problem  were  described  in  Section  1 of 
Chapter  III  and  they  are  summarized  here. 

1.  The  signal  Is  a dc  signal  with  value  10.0. 

2.  The  noise  was  generated  using  the  computer's  uniformly  dlstribu 
ted  random  number  generator  and  transforming  to  zero  mean  Gaus- 
sian random  numbers  with  the  Box-Muller  equations. 

3.  The  initial  weight  vector,  W(l),  for  the  LMS  algorithm  was  set 
equal  to  5.0. 

4.  The  mean  square  error  was  computed  using  an  iterative  estimate: 


(89) 
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• si>2 

ECtij  - »/]  = i^-j (o.i) 

5.  The  mean  square  error  Is  computed  for  500  Iterations  to  produce 
one  ensemble  curve.  Twenty-five  ensembles  are  averaged  point 
by  point  and  this  curve  Is  plotted  as  a function  of  the  number 
of  iterations. 

Program  Structure: 

A separate  program  Is  employed  for  each  algorithm.  The  program 
ESTSIG  implements  the  expected  value  algorithm  and  the  program  ANEST  imple- 
ments the  LMS  algorithm.  Each  program  is  divided  Into  two  main  parts.  One 
part  generates  the  noise,  computes  the  weight  vector,  and  determines  the 
mean  square  error.  The  other  part  averages  the  ensembles  and  plots  the  re- 
sult. Figures  0.1  and  0.2  show  the  structure  of  the  programs. 

■roqram  Segment  Descriptions 

1.  ESTSIG  — Main  program  to  implement  the  expected  value  algorithm  and 

compute  the  mean  square  error;  calls  subroutines  CGAUSS 
and  PLOTIT. 

INPUTS:  RO  — correlation  between  the  two  noises 

2 

XVAR  — variance  of  the  noise  n^  E[n^  ] 

2 

YVAR  --  variance  of  the  noise  n2;  Efr^  ] 

OUTPUTS:  plot  of  the  mean  square  error  as  a function  of  the  number 

of  Iterations. 

2.  ANEST  — Main  program  to  Implement  the  LMS  algorithm  and  compute 


the  mean  square  error;  calls  subroutines  CGAUSS  and  PLOTIT. 
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INPUTS:  C — term  which  determines  the  rate  of  convergence- 

RO  — correlation  between  the  two  noises. 

2 

XVAR  — variance  of  the  noise  n^ ; E[n^  ]. 

2 

YVAR  — variance  of  the  noise  n^’,  ECn^  3 . 

OUTPUTS:  plot  of  the  mean  square  error  as  a function  of  the  number 
of  Iterations. 

3.  CGAUSS  (II ,I2,R0,X,Y)  — generates  a pair  of  zero  mean  Gaussian  ran- 
dom numbers  using  the  Box-Muller  equations 
to  transform  two  uniformly  distributed 
random  numbers;  calls  library  uniformly 
distributed  random  number  generator  RANDU. 

— seed  number 

— desired  correlation  between  the  pair  of  numbers. 

2 

— desired  variance  E[n^  ] through  COMMON  statement. 

2 

— desired  variance  Efr^  ] through  COMMON  statement. 

— zero  mean  Gaussian  random  number 

— zero  mean  Gaussian  random  number 

— seed  number. 

4.  PLOTIT  (X,NP,NV,NS)  — plots  mean  square  error  as  a function  of  the 

number  of  iterations;  calls  library  plotting 
subroutine  TEKPIT. 

INPUTS:  X — matrix  where  each  row  Is  a vector  of  points  to  be 
plotted  so  that  up  to  five  curves  can  be  plotted 
on  the  same  graph;  can  also  be  entered  as  a vector 
If  only  one  curve  Is  to  be  plotted. 


INPUTS:  II 
RO 

XVAR 

YVAR 

OUTPUTS:  X 
Y 
12 


NP  ~ number  of  points  within  each  vector  for  plotting- 
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NV  — number  of  variables  or  curves  to  be  plotted  on  the 
same  graph. 

NS  — scaling  determinant;  value  of  zero  allows  programmer 
to  scale  while  value  of  one  provides  automatic 
scaling. 

OUTPUTS:  plot  of  vectors  entered  as  a function  of  the  number  of 
points  within  the  vector. 


B.  Three  Inout  Noise  Canceling  System 


This  computer  simulation  extends  the  two  input  system  in  order  to 
further  compare  the  two  algorithms.  The  equations  for  the  two  algorithms 
are  given  in  Figures  4.4  and  4,5.  These  were  implemented  on  the  POP-11/45 

I 

computer  at  Duke  University.  The  mean  square  error  was  plotted  as  a function 
of  the  number  of  iterations. 


Specifications: 

The  specifications  for  the  problem  were  described  in  Section  1 of 
Chapter  V and  are  summarized  here. 

1.  The  signal  Is  a dc  signal  with  value  10.0. 

2.  The  noise  was  generated  with  a method  developed  by  Charles  Liu 
and  described  In  Appendix  B.  Independent  Gaussian  random  numbers 
are  transformed  to  correlated  ones  using  this  method. 

3.  The  initial  weights  for  the  LMS  algorithm  were  set  to  arbitrary 


values:  W^( 1 ) ■ 5.0  and  W3(l)  * 3.0. 


4.  The  mean  square  error  was  computed  using  an  iterative  estimate: 
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ESTSIG 

C IST3IG  **  TO  ESTIMATE  A SIGNAL 

C 

-DIMENSION  SIG(300), TEMP (800), TMSEl 300), TM3E2 (500) 
COMMON  /A/XVAR#YVAR 
'DATA  11,12/16401,0/ 

KRITE  (6,90) 

READ  (6,91)  RQ 
WRITE  (6,30) 

REAC  (6,83)  AVAR 
NRIT6  (4,81) 

:REA0  (6,33)  YVAR 

iN0ITER*2S 

:NS«0 

•3 IGNAL*  10 ,0 

!NR ITE  (6,94)  RO,  XVAR , YVAR 
'CO  iee8  JKa  1 ,300 
1000  TMSE ( JK)«0.0 

■DC  23  1*1,300 
29  S1G(I)*3IGNAL 

•DC  180  K*  1 , rNOZTER 

LI *0,0 

L2*0,0 

CC  101  JKL*  1 , 900 
101  T^SE2 ( JKL)*Q«9 

CC  10  1*1 , Sv)C 

CALL  CGAUSS  ( I l , 12, RC , X , Y ) 

Xl«9S6<n  + X 

X2*Y 

Lt*UUXl«X2 

•L2*U2*X2*X2 

13*U1/U2 

TEMP(I)*((X1«U3*X2)-S1G{I))**2 
CO  11  J*1,I 

11  Tm3E2(I)*T:M3E2(I)  ♦ TEMP  < J) 

TM3£2(I)*TMSE2(I)/(PLCATCI)) 

19  TMSE(I)*TM3E(I)  ♦ TM362 ( 2 ) 

100  CONTINUE 

CC  2e0  1L*1,S00 

200  TM3E(IL)*TM3E(IL)/(fLCAT(N0ITER)) 

CALL  PLOTIT  (TM3E, 322,1, N3) 

C 

00  FORMAT  (/IX,  »XVAR»,/1X,  ♦ 9«509cccccfl«<j?ccc  * ) 

01  FORMAT  (/IX, *YYAR»,/1X, • »C699C9C990«tf •• » ) 

03  FORMAT  (cl«,S) 

04  FORMAT  (/iX,»NOIT£R»,/iX, t «««««!) 

00  FORMAT  (15) 

90  FORMAT  (/IX,  »RC»,/1X,  • oriocsptcoe ») 

91  FORMA?  (FIS ,3) 

92  FORMAT  (/lX,<N3»,/iXf »e«) 

93  FORMAT  (ID 
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94  FORMAT  <lX,P10t5,lX,ei«t8,iX,E!6.8) 

9fl  FORMAT  C/U/ ’SIGNAL*, /IX,  » oo«oo<»oo«»  j ) 

97  FORMAT  <FU,n 

3008  CAUL  -EXIT 

■ I NO 


id 


ANEST 
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ANEST  IS  A PROGRAM  TC  ESTIMATE  A SIGNAL 
USING  THE  .LN3  ALGORITHM 

DIMENSION  M501),  S IG (30e ) , TMSE  (300)/  TMSE2 (500),  TEMP (500) 
COMMON  /A/XYAR,YVAR 
CATA  11,12/16401,0/ 

WRITE  (6,93) 

READ  (6,91)  C 
WRITE  (6,93) 

READ  (6,83)  RQ 
HRITE  (6,85) 

READ  (6,91)  XYAR 
WRITE  (6,86) 

READ  (6,91)  YVAR 

:NQITER«1 

■N3»0 

S ISN  AL*  1 0 ,8  1 

WRITE  (6,80)  C,  RO,  XVAR, YVAR 

CO  1000  L*  1 ,330 

TM3E(L)*0t8 

CC  25  1*1,303 

SI6(I)«SIGNAL 

CO  1 20  K*i ,NOlT£R 

CO  121  JKL* 1 , 500 

TMSE2 ( JKL) *3,0 

W(JML}«0,0 

W ( 1 ) »5 , 3 

'CO  68  1*1,500 

CALL  CGAUSS  (II, I2,P',X,Y) 

Xl*31GU)*X 

X2«Y 

.3MAT«XlfWC)*X2 

W(I*i)*M(l)42,a*C*X2*C3HAT) 

TEMP(I)«(3MAT»3lG(I))*«a 
CO  11  J*t,I 

TMSE2(1)»TM3E2(I)  ♦ TEMPfJ) 

TMSE2<I)«TM2f2(I)/C?LCAT(I)) 

TM3E (IMTM8EO)  ♦ TMSE2U) 

CONTINUE 
CO  200  ML*  t ,530 

TM3E (KL) *TM3E (KL) /(FLOAT (NCITER)) 

CALL  PLOTZT  (TMSE, 588,1, NS) 

FORMA" (13) 

•FORMATt/lXf  »NOmR,#/lX»,«f9*9l  > 

FORMAT  (/IX, Ell, 8, 1?, F 13,5, IX, 2 (El 6, 8, IX)) 
format </ix » » signal', /ix, »«9««f «««««»} 
fcrmat(f;«,i) 

FORMAT (F10j6) 
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88  FGRFAT  C/1X  * »XVARI,/IX,  » I) 

0<S  FORMAT  (/lKr » Y V AR  * ,/lX,  • ocoeo«e9«eoco*«fl  » J 

98  FORK AT (/l X# ' X C l ) « #/lX, «««« «««««« 

91  F0RMT(Ei6,8) 

92  FCRHAT(/1X,  » C » # ✓ l X # » COq<»oo«ot:«o«C*«« » ) 

93  FQRMAT</tX, ‘ROt,/lX, * 

96  FCRFAT</lX,»N3»,/iX, »«») 

97  FCRHATUO 

3000  CALI.  EXIT 

ENO 
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SUBROUTINE  C0AU3S  (11/ 12/ FO, X, Y> 

TO  GENERATE  CORRELATED  GALSSIAN  RAND CP  NUMBERS 

CIMEN3I0N  A (2) 

COMMON  /A/X7AR/YVAR 

CATA  PI, XMU/TMU/^, 14155, 0,,0,/ 

XXVAR«8SRT(XVAR) 

YYVAR*5GR7 ( YYAR) 

CO  20  J«l,2 

CALL  RANOU  (Z1/I2/U) 

IF  (U  ,EQ,  1,00000)  CC  TO  10 
IF  (U  ,£Q,  0,00000)  CC  TO  10 
A ( J) *U 
CONTINUE 

C*S0RT(«2»3*ALQG(A(1))) 

C«3IN(2,0*FI*A(2)) 

E«(SGR7(1,C-RO*RC))MCC3(2,0*PI*A(2)))  + RC*C 

X*XMU+XXVAR*C*8 

Y»YMU+YYYAR*C*C 

RETURN 

eno 
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P10TIT 


SUBROUTINE  PlCTIT<X,NF,NY,N5) 

DIMENSION  X^V>NP),XMAx<S),mM3),XJ'EAM5) 

SCALE  CATA 

DC  30  ’J*t,NY 
XFEANCJ)W0,0 
XHAX(J3*0t9 
CO  10  1*1, NP 

X«ANCJ)*XNEAN(4)tX(4,n/fLCAT(,NP) 

10  XNAX(J)*AI*aXI  tXNAX(  J3fXCJr  X)) 
XMN(J3*X*AXf4) 

DC  IS  1*1, NP 

13  XPIN<J)*AMIN1(XNIN(J),XCJ#I)) 

IF  CXHAX ( J) } 22,22,2 3 
2e  XMX{J)»XMIMJJ 
DO  23  1*1, NP 

23  XKAX(J)*A«.4X1(XKAXC.J),X(J,D) 

30  CONTINUE 

XFAXX*XMAX(13 

XPINI«XNIN<1) 

CO  27  1*1, NV 

XHAXI»AMAXHXNAXI,XMAX{I)) 

27  XPIKX«AHlNUXPlNIrXMlNCI)) 

Z?(S3*£G*0)  GG  TO  30 
YlNT»(Xi*AXI*»>yINI7/2a4>FINI 
40  IF(ABS(XMAxI*TlNT) #NS,3v2)  SO  TC  4« 

YS§1, 

80  TO  50 

43  rS*313,/ABSU?(AXJ®TlNT3 
50  X3*300,/NP 

Ir(NS,Sc,3)  SO  TC  S3 
CALL  TEKPLT13, 0,0,3, 03 
WRITE (5,2322) 

* DC  59  1*1, NV 

«0  *SJU<6,2400)  X»XNAXU)»XNlMI),XKEAhU) 
WRITE (5 , 2360  3 
REACU,  1 1303  Y3 

!F1ABS<YS*XNAXI3  ,6Y,32,E44)  YS*32,S4<5m AX  I 

YS*4fl,/TS 

MR  Z ?E ( 6 ; 2000} 

REAC(6,119S)  YINT 
03  CALL  TcNPLTlS, 0,0, 0,33 

CALL  TEKPLT(2,30,10>3S3,1S) 

CALL  TEKFLT12, 33, 5^3,550, 540 3 
CALL  TEKPir(2,30,if*,3e,$40J 
CALL  TEKPLT12,830,10,S*e,6433 
CALL  TEKFL7C, 55, 323,850,3253 
DC  100  J»1,1S 
1*949, «(I3«J}*10, 
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CALL  TEKPLTC2#S0,  IY»50,XT> 

CALL  TEKPLTt2»840flYf850fIY) 

IPfYINTf<45,/Y3)MLCA7<«*J) 

IYiIY*6 

CALL  TEKPLU3#85SfIY>ef8) 

108  MPITE<8,2000)  X? 

00  190  J»l*18 
1X«90*J490 

CALL  TEKPLT(2r IX,18,IX,20) 

CALL  TfKPLTt2fXX#320/XX/330) 

158  CALL  TEKPLT(2fZXf630flX?<40) 

0C  600  J'1,*V 
YT*<X(J,1)*YINT)M3*325# 

XP<A83(YY),Ce,3flE*6)  YY*3«£46MY/AB3(YY) 

I YY»XNT(YY) 

CALL  TEKPLT<3#5,XYY#5, JYY) 

HfilTE<5,2700)  3 
00  808  I«2,NP 

U*(X(J,  J-U*YIMT)«YSt325, 

IF(Ae3(Yl)#6E#3,l£*8)  Y1*3«E+€*Y1/AB3(Y1) 

IYl«IhT(Yl) 

Y2«CX<0,I)-YIHT)*Y3*355, 

XP(A£8(Y2)t6Et3tlE*«)  Y2«3 , E+8* Y2/A88 (Y2) 

IY2»InT(Y2) 

IX1»(I*XS)+B0, 
lX2«(!tiJ*X3.»9e# 

808  CALL  TEKPLT12, XXI , I Yt> 1X2, IY2) 

1006  FCRMAT<1X,I8) 

1106  PCRMAT(iX,Fl«,0) 

2000  FCf»MTUH*,lFE10,3) 

2300  FORMAT (/IX,  «VAPIAeiE»f5X,»YFAX»#nXf  » YUM  , liX,  » MEAN’//) 
2466  FCRMAT(lX,l3,3(5X,l?Ei€«3)) 

2500  FORMAT (//1X«  ’ENTER  Y-AXI3  SCALE  » /1XI?  » < )’) 

2000  FGRMATC//IXM  ENTER  X*AXJS  PC3I7ICM/1X,  * t )») 

2700  FORMAT ( lH4f |2)  ‘ 

2200  FORMAT (///) 

RETURN 

ENC 


H<s‘i  - s<>2 

EC(Jj  - Sj)a3  s — 3 (0.1) 

5.  The  mean  square  error  Is  computed  for  500  Iterations  to  produce 
one  ensemble  curve.  Twenty-five  ensembles  are  averaged  point  by 
point  and  this  curve  is  plotted  as  a function  of  the  number  of 
iterati ons. 

Program  Structure 

A separate  program  is  utilized  for  each  algorithm.  The  programs 
ESTSI2  and  ANEST2  implement  the  expected  value  and  IMS  algorithms,  respect- 
ively. Each  program  is  divided  into  three  main  parts.  The  first  part  finds 
the  eigenvalues  and  eigenvectors  of  the  input  correlation  matrix  so  that 
correlated  Gaussian  random  variables  can  be  generated.  The  second  section 
generates  the  noise,  computes  the  weight  vector,  and  calculates  the  mean 
square  error.  The  third  section  averages  the  ensemble  curves  and  plots  the 
resulting  curve.  Figures  0.3  and  0.4  show  the  structure  of  the  programs. 

Program  Segment  Descriptions 

1.  ESTSI2  — Main  program  to  implement  the  expected  value  algorithm  and 

compute  the  mean  square  error;  calls  subroutines  JACOBI, 
GAUSS3,  and  PIOTIT. 

INPUTS:  VARNI  — variance  of  the  noise  n^ ; E[n.^]. 

OUTPUTS:  plot  of  the  mean  square  error  as  a function  of  the  number 

of  Iterations. 

2.  ANEST2  — Main  program  to  implement  the  IMS  algorithm  and  compute 

the  mean  square  error;  calls  subroutines  JACOBI,  GAUSS3, 
and  PLOT IT. 
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INPUTS:  VARNI  — variance  of  the  noise  n-j ; E[n^]. 

CMU  — rate  of  convergence  term,  u. 

OUTPUTS:  plot  of  the  mean  square  error  as  a function  of  the  number 

of  iterations. 

3.  JACOBI (T,  EIGEN,  R012,  R013,  R023) 

— subroutine  to  compute  the  eigenvalues  and  eigen- 
vectors of  a symmetric  matrix.  This  program  is 
a variation  of  the  one  listed  in  Applied  Numeri- 
cal Methods.  It  is  used  to  transform  Independent 
Gaussian  random  numbers  to  correlated  ones. 

INPUTS:  A(l,2)  — correlation,  between  the  noises  n-j  and  n 

A(l,3)  — correlation,  p^»  between  the  noises  n^  and  n^. 

A(2,3)  — correlation,  p23,  between  the  noises  n2  and  n^. 

OUTPUTS:  T — matrix  where  each  column  is  an  eigenvector. 

EIGEN  — vector  of  eigenvalues 

R012  — correlation,  p^»  between  the  noises  n^  and  n 2* 

R013  — correlation,  p^»  between  the  noises  n-|  and  n3« 

R023  — correlation,  p23,  between  the  noises  n2  and  n3> 

list  of  the  elements  of  the  input  matrix,  A. 

4.  GAUSS3  (IX,S,AM,F)  — subroutine  to  generate  three  Independent  Gaus- 

sian random  numbers;  calls  GAUSS. 

INPUTS:  IX  — seed  number 

S — standard  deviation  of  the  variance  of  each  random 
number;  7^  • 

AM  — desired  mean  of  the  generated  Gaussian  random 
numbers. 
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OUTPUTS:  F — vector  containing  three  Independent  Gaussian 

random  numbers. 

5.  GAUSS  (IX,  S,  AM,  V)  — scientific  subroutine  package  program  to 

generate  one  independent  Gaussian  random 
number;  calls  library  subroutine  RANDU. 

INPUTS:  IX  — seed  number 

S — standard  deviation  of  the  desired  random  number. 

AM  — desired  mean  of  the  random  number. 

OUTPUTS:  V — Gaussian  random  number  generated. 

6.  PLOTIT  (X,  NP,  NV,  NS)  — plots  mean  squares  error  as  a function  of 

the  number  of  iterations;  calls  library 
plotting  subroutine  TEKPLT. 

INPUTS:  X — matrix  where  each  row  is  a vector  of  points  to 

be  plotted  so  that  up  to  five  curves  can  be 

j 

plotted  on  the  same  graph;  can  also  be  entered 
as  a vector  if  only  one  curve  is  to  be  plotted. 

I 

NP  — number  of  points  within  each  vector  to  be  plotted. 

NV  — number  of  variables  or  curves  to  be  plotted  on 

the  same  graph. 

NS  — scaling  determinant;  value  of  zero  allows  pro- 
grammer to  scale  while  value  of  one  provides 
automatic  scaling. 

OUTPUTS:  plot  of  vectors  entered  as  function  of  the  number  of  points 
within  the  vector. 
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:EST3I2  — TO  ESTIMATE  A SIGNAL  USING  THE  EXPECTED 
VALUE  ALGORITHM  ANC  TNG  NOISE  'ONLY 
INPUTS 

DIMENSION  3ieC60e>,7!MPC5e0),TMSE(30fi),TNSE2(500) 
DIMENSION  CIAG(3),Y(3)#VECS(3,3),VALS(3), F (3) , EDEN (3,3) 
DIMENSION  EI3EN(3)#T(3,3),TRANS(3,3) 

DOUBLE  PRECISION  VEC3, V ALS , DSC  12, 0RC1 3 , OR023 
DATA  IX, AH/16401, 0t£/ 

NRITE  (4,80) 

REAO  (6,83)  VARN1 
NRITE  (6,81) 

8*3CRT (YAfiM  ) 

N0ITER»25 
•DC  1000  JM 1 ,900 
•318(JK)>10(0 
TM8E ( JK) *0 v0 
-DC  10  1*1,3 
CO  10  J*l,3 
EDEN ( I , J) *0 ,0 

CALL  JACOBI  (VEC3,VALS,GRQ12,CRC13,CRC23) 

DC  22  1*1,3 

EI6ENCn*3NGL(VAL8(J)) 

CO  20  J*  1 ,3 

TCI,J)*SNGL(VEC8(!,J)) 

RC12*SNGL(CRC12) 

FC13*8NGl(CRC13) 

R023*SNGL(CRC23) 

CO  21  1*1,3 

C1AG(1)«8SRTCEXGENU)) 

DC  22  1*1,3 
EDEN(I, 1)*0IAG(I) 

CO  IS  1*1,3 
-CO  18  J*1 , 3 

1RAN8(I,J)*T(I,1)*ECEN(1, J)AT(I,2)*ECENC2, J)+ 
T(I,3)«E0EN(3,J) 

CO  100  K*1,MC1TER 
A I2*f  , 0 

M3>M 

A22*«,00t 

423*9,0 

A33*0,8 

CO  ft  I JKLf 1,300 
TMSE2 ( JKL) *0«0 
DO  2C2  1*1,600 
CALL  GAU3S3  (IX,S,AM,P) 

DO  17  N* l , 3 

Y(N)*TRAN$(N,i>*F(l)<TRAN3(N,2>*F(2)-»TRANS(N,3>*F(3) 

X 1 *S!G ( 1 )* Y ( 1 ) 

X2*Y (2) 
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X3*Y<3> 

A12*Al2>Xl*X2 
At3«Al3m*X3 
A22*A22+X2*X2 
A33?A33+X3*X3 
.A23*  A23*X2*X3 
kOEN* A22*A33"A23* A23 
F2«(Al2*A33-Al3*A23mCEN 
»3«(A13*A22-A12*A23)/KCEN 
TENP(J)«(CX1**2*X2«I*:*X3)-SIG(I))**2 
CO  it  J«1  # ! 

11  TM3E2 ( 1 ) sTy  SE2  ( I ) ♦ TEFPtJ) 

TP3E2(I)»T)»SE2(I)/(FICAT(I)) 

2B0  TNSE(I)«TFSECI)  t Tl'JE2<I) 

;aa  continue 

CO  291  11*1/999 

291  THSECIL)«T*3ECIL)/CHCATCNCXTER)) 

CALL  PLOT IT  ( TNSE# 320 , 1 , 0 ) 

C 

89  FCRNAT  (/IX  / 1 V ARM  *// IX* » 0009Cff«9«  » ) 

31  FORMAT C/8X, »A<1,2) »/ ' A ( l , 3 J » , U (2/3) t ,/6X, » OC? 39* » , 

1 » 9*3933 I , I 993*99  » ) 

83  FCRKA7  CE18,3)  ' 

90  FORMAT  { F3 , 1 ) 

3009  CALL  EXIT 

ENO 
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ANEST2 

ANE3T2  *-*  TO  ESTIMATE  A SIGNAL  USING  THE  LM3 
ALGORITHM  and  two  ncise  only  inputs 


•DIMENSION  JIG(300),TEMP(5B0),7MSE(50(!),7M2E2(500) 
DIMENSION  0IAG(3),Y(3),VEC3(3,3),VALS(3),F(3),EDEN(3,3) 
DIMENSION  EIGEN  (3 ) #7 (3 , 3), TRAN S (3, 3 ) , N2 CEB  1 ) , H3  (301 ) 
DOUBLE  PRECISION  VEC 8 , Y AL3 , DSC  1 2 , CR013 , CF023 
DATA  IX, AM/18401, 3te/ 

HRITE  (6,80) 

READ  (8,83)  YARN) 

WRITE  (8,84) 

READ  (6,83)  CMV 
NRITE  (6,81) 

3-3GRT  (VARM) 

NOJTSF-25 
DC  1000  JK« 1,800 
.SICCJK)*10,0 
1000  TM3E ( JK) -0*3 

CO  10  1-1,3 
DC  10  J = l,3 
18  EDEN ( I , j) afi  ,0 

CALL  JACQ3I  5VECSrVALS,CRC12,CR013,0RC23) 

DC  28  1-1,3 

e:gen(Z)-3ngl(valsu)) 

DO  22  0*1,3 

20  7(I,J)«SNGL(VEC3(I,;)) 

R012«8N5L(CR012) 

R013«3N€L(DR013) 

RQ23«3NGL(CRC23) 

DO  21  1-1,3 

21  DIAG(I)-36RT(EIGEN(1)) 

CO  22  1-1,3 

22  EOEN(I,  I)-CIAG(I) 

DC  IS  1-1,3 

00  13  J-1,3 

15  TRANS(l,J)*TCI,i)*EDEN(l,J)M(I,2)*ECEN(2,J)H 

* T ( 1,3) -EDEN (3, J) 

•DO  100  K- 1 , NO  ITER 
DO  101  JNL*1,608 
N2( JRL)-0»0 
N3(JNL)-0.e 

101  TM3E2 ( JKL)-0«0 

*2(l)-3,0 
N3 ( 1 ) -3,0 
DO  200  1-1,900 
CALL  CAUS33  (IX, 8, AM, P) 

DO  17  N-1,3 

Y(N)-TRAnS(N,1)*F(1)1TRANS(N,2)-P(2)yTRAN3(N,3)-F(3) 
.Xl-SIG  (I)TY  ( 1) 
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X2«Y(2) 

>3-Y  <3  J 

■ 3HAT«X1»W2(I)*X2«>H3(I)*X3 
X2<IM>«*2(I)*2,0*CPUa3HATAX2 
*3(IM)**3<J)+240*CK*3HA7*X3 
TEMP(I)«(3VAT-SIG(I))**2 
•CC  11  J»l/I 

11  TP3E2(I)*TM3E2(I)  ♦ TEPPCJ) 

TP3E2(I)«TP3E2(l)/(FlGAT(I)) 

289  TK5E(I)«TM3E(I)  ♦ XPCS2(X} 

138  CONTINUE 

■CO  221  IL  » 1 ; 339 

291  TM3E(IL)»TME<Il)/CFUA7(S0I7ER)) 

• CALL  PIOTIT  <TM3E,52g,i,0) 

C 

88  FORMAT  (/IX,  ' VARNt » , / tX,  » 9999*  9*9990*9999 » ) 

81  FORMAT (/fiX* »• A ( 1 , 2 ) * , I A ( 1 , i ) * , » A (2, 3 } l , /€X , * 99**9*  » , 

1 »99909*  » , » 99*9991 ) 

83  FORMAT  (El«,8) 

84  FORMAT  (/IX, * CPU  * , /IX  , » 90e*9**0***99e9C » ) 

90  FORMAT  (F3,n 

3938  CALU  EXIT 

END 


•rvr»-o>  cx  - n o m >•»  o o ooo 
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JACOBI 


SUBROUTINE  JACOBI  CT , E ! 6EK / RC 12, RC13, RC23 ) 
EIGENVALUES  ANO  EIGENVECTORS  EY  TNE  JACCBI  METHOD 
APPLIED  NUeERICAL  HETHCOS,  EXAMPLE  4,4 

IMPLICIT  REALMS  (A*M/(0«Z) 

DIMENSION  A <3, 3 ),T (3,3), A JK (3). EIGEN (3) 

DATA  N» ITMAX/3* 10/ 

DATA  EPSi, EP32#EpS3/l,00C!-ia,l,00C-i0  / 1,280^5/ 

READ  PARAMETERS  AND  :EST.ABU3M  STARTING  MATRIX  A 
CO  2 I«I,3 
•CO  2 J«l,3 
T(i/J)«a,e 
A(I,J)«0t8 

READ  (0/ 182}  A{i,2)/A<l/3)/A(2/3) 

A(2/l)aAU#2) 

AC3,l)«A(l,3) 

A(3,2)*A(2,3) 

ACt/t>»l098 
A(2,2)«l,08 
A (3/3) 3l (88 
RG12*A (1/2) 

FCi3«A(i/3) 

RG2J« A (2/3) 

N.RITE  (fi/191)  ((ACI/J)rJ<«l.»N)/IflrN) 

NM1=«M 

SET  UP  INITIAL  MATRIX  T/  COMPUTE  SIGMA!  AND  3 
SIGMA!  *0,3 
iCFFCSG«0,3 
CO  3 I * I / N 

•SIGMAl*8JGMAlAA(J/I)AA3 
T (1/ 1) *1,3 
IPJ*U1 

IF  (X  ,GE,  N)  GO  TO  0 
CC  5 «« IP  1 / N 
CFFDSC*CFPC8fl'»A(I,J)«A2 
S*2,0*CFFO8GASIGMAI 

•BEGIN  JACOBI  ITERATION 
DC  20  ITSR«1/ITMAX 
■CO  28  IM/NMi 
IP  1 *•  1 ♦ 1 

-CC  28  J«TPI,N 
CSDAES(A(J/I)*A(J/J)) 

C COMPUTE  SINS  ANO  COSINE  OF  ROTATION  ANGLE 

IF  (C  ,LE«  epsu  GO  TO  9 
IF  iCABS(A(I/4J)  »LE * EPS2)  SC  ?C  20 
F42,8»a(1/J)«q/ou/1>*a<j/J)) 
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.3PQ*C3CRT(F*PtQ*«) 

CSA*C3GRTUl*0*fl/3P6)/2*9) 

.3nA*P/<2,0*C3A*SP«) 

•60  TC  10 

9 C3A«ita/o3GRTC2taoa; 

3NA*C3A 

19  CONTINUE 

C ' 

0 LP0ATE  C0LUH3N  I ANC  J CP  T * ECUIYAUM  TO 

C -MULTIPLICATION  BY  Xh 5 ANNIHILATION  MATRIX 

-CO  11  K«1,N 

tMQLCXI*T(K,  I) 

T(X#I)*NQLCXI*C3A*T(X/«1)*3NA 
11  7(K,J)*HOLCkI*3HA-T(X,-J)*C3A 

C 

C COMPUTE  NE*  ELEMENTS  >CP  A IN  ROMS  I ANC  J 

CO  IB  X * I / N 

IP  <X  *GT*  J)  60  TO  13 
AIK<X)*A<X,K) 

A (X, X) »C3A* A XX (X)+3N**A(X,J) 

IP  (X  *NE*  ‘J)  60  TO  U 

A(J.K)*3KAAAZ:«<X)*C3A*A(J,XJ 
H 60  TO  14 

15  !NCLCXK*A(I#X) 

, ACI^K)«C3A*HCLCIX+3NA*A<J#X) 

A(.J,K)»3N.AAHQL2XX*Ci;»A  (J,X) 

14  -CQNflNUS 

C' 

C COMPUTE  NE*  ELEMENTS  CF  A IN  COLUMNS  I AND  J 

aikcj>*snA!MIkcx)i.csa*ajx;J) 

c 

C XHEN  X 13  .LARGER  THAN  I 

-CO  19  X*l,\l 
IF  (X  «LE«  I)  GO  TO  14 
•A(K,  J)*3NA«AIXCX5-C3aAACX,J) 

•GO  TO  19 

13  hOLCKI«A(X#n 

A(X,  I)«C3AapCLCKI*3M*A<Kf  JJ 
-A  (X«  J)  *3NA*rtCLCK I«C3 A * A (X , J) 

19  CONTINUE 

20  A<I/\J)*fl,8 
C 

C PINO  3IGMA2  POP  TRANSFORMED  A ANC  TEST  FOR  CONVERGENCE 

.3IG*A2*0t0 
CC  21  1*1 rN 
EIGEN ( Z)*A( I » 1) 

31  3IGMA2*3!CM.A2+E!SEN\3)a*2 

IF  Cl*0-3IGXAl/8IGHA2  ,rGE.  E333)  GO  TC  .24 
GO  TO  300 

26  -SIGMA  \33XGMA2 

C 

c IP  ITER  EXCEEDS  XTNAX#  NG  CONVERGENCE 

* R I T E (6,223)  ITSRiS^IGMAUSlGNAa 


»-*•  0"0 


01 

102 

203 


C 

300 


* 

-* 

:* 


EROM  CVtY  ffURftiiSHF.fi 


40  TO  1 

FORMATS  FOR  INPUT  ANC  .'CUPtT  •.STATEMENTS 
FORMAT  (18X,«F10.4) 

FORMAT  (5X,3F4,4) 

FORMAT </5X#  • NO  CONVERGENCE#  UTH'»//7X# 

* ITER  * M3,3X#»  3 :«  *,  FIB. 5/ S>, 

» 3I5MA1  • »#Fl8,3fSX,*  SI6MA2  • » , F 10 *5 , //SX , 
» THE  CURRENT  A MATRIX  IS*#/!*) 

RETURN 

END 
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GAUSS3 


SUBROUTINE  GAU3S3  (UrS,AN,R) 

TO  GENERATE  3 XNOERLENCENT  GAIS3IAN  RANDOM  NUMBERS 

•DIMENSION  F(3) 

•CC  12  I*l,3 
CALL  GAU03  (IX, 3, AM, Y) 

19  E ( I) *Y 

RETURN 

:ino 
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12.  All  TRACT 


The  purpose  of  this  study  is  to  compare  two  noise  canceling  system  algorithms  when 
the  system  has  a small  number  of  inputs.  The  expected  value  and  LMS  (least  mean  square) 
algorithms  attempt  through  different  methods  to  determine  the  optimum  solution  given  by 
the  Wlener-Hopf  weight  vector.  The  expected  value  algorithm  estimates  the  necessary  ex- 
pected values,  substitutes  them  into  the  Wiener-Hopf  equation,  and  performs  the  necessary 
matrix  Inversion.  The  LMS  algorithm  uses  gradient  techniques  and  the  method  of  steepest 
descent  to  avoid  the  complexity  of  the. matrix  inversion.  The  algorithms  are  compared  for 
two  and  three  input  systems.  In  terms  of  the  rate  of  convergence  of  the  mean  square  error 
and  the  minimum  error  after  five  hundred  iterations. 

A computer  simulation  of  a dc  signal  in  Gaussian  noise  provided  the  comparison  be- 
tween the  two  algorithms.  The  LMS  algorithm  presented  problems  in  implementation  because 
of  two  parameter  values  that  had  to  be  chosen.  Although  a value  for  the  initial  weights 
had  to  be  determined,  this  value  was  not  as  important  as  was  the  value  of  u,  the  term  that 
controls  the  rate  of  convergence.  This  term  must  be  within  a specific  range. in  order  to 
guarantee  convergence  of  the  algorithm. 

In. conclusion,  the  expected  value  algorithm  was  judged  to  be  the  better  of  the  two. 
The  tradeoff  in  terms  of  performance  for  simplicity  In  the  LMS  algorithm. was  great  es- 
pecially when  the  slgnal-to-nolse  ra.t1o  or  the  correlation  between  the  noises  was  small.  • 
If  a "correct"  value  of  u Is  chosen,  the  results  may  be  similar,  but  if  not,  the  results 
can  be  extremely  different  from  those  of  the  expected  value  algorithm.  Thus  for  a small 
number  of  Inputs,  where  the  matrix  to  be  Inverted  is  of  small  dimension,  it  is  bettec  in  . 
terms  of  performance  to  use  the  approximation  to  the  optimum  solution  given  by  the  ex- 
Ipected  value  algorithm. 
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